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Abstract 


This  thesis  examined  the  suitability  of  a  particle  simulator  known  as 
“Monte  Carlo  n-Part-icle  6  (MCNP6)  to  determine  how  well  U.S.  Army  vehicles 
shield  their  occupants  from  gamma  radiation.  MCNP6  compared  favorably  to 
three  alternate  modeling  approaches  and  another  particle  simulator.  In  a 
validation  experiment,  MCNP6  was  found  to  produce  a  gamma  protection  factor 
within  5%  of  the  experimental  result  with  95%  confidence.  Further  evaluations 
validated  MCNPG’s  ability  to  produce  reliable  results  with  simplified  or 
approximated  inputs,  model  different  sources,  and  position  those  sources  in 
different  shielding  geometries.  MCNP6  is  suitable  for  evaluating  gamma 
protection  of  Army  vehicles  and  further  research  to  explore  the  effects  of  different 
shielding  shapes  and  evaluate  the  suitability  for  MCNP6  to  calculate  a 
comprehensive  radiation  protection  factor  is  warranted. 
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VERIFICATION  AND  VALIDATION  OF  MONTE-CARLO  N-PARTICLE  6  FOR 


COMPUTING  GAMMA  PROTECTION  FACTORS 

1.  Introduction 

1.1.  Motivation 

How  well  does  the  United  States  Anny’s  current  fleet  of  combat  vehicles  protect 
soldiers  from  the  unique  effects  of  nuclear  weapons?  This  simple  question  underscores  a 
wide  variety  of  tactical,  operational,  and  strategic  concerns.  Exploding  nuclear  weapons 
emit  powerful,  invisible  pulses  of  ionizing  radiation  in  addition  to  their  devastating 
thennal  and  blast  effects.  Unlike  the  blast  and  thennal  effects,  the  ionizing  radiation 
hazards  persist  long  after  the  explosion  in  the  form  of  radioactive  fallout.  Both  the  initial 
radiation  pulse  and  the  fallout  can  cause  radiation  sickness  and  death  in  the  short  tenn.  In 
the  long  term,  heavy  doses  of  ionizing  radiation  can  cause  pennanent  health  problems 
and  many  forms  of  cancer. 

In  order  to  manage  these  risks,  the  Anny  needs  reliable  and  verifiable  information 
on  how  well  their  vehicles  protect  crews  and  passengers  from  ionizing  radiation.  For 
planning  purposes,  this  infonnation  is  expressed  as  a  radiation  protection  factor  (RPF) 
that  compares  the  expected  radiation  effects  on  an  unshielded  subject  to  the  expected 
effects  on  a  shielded  subject  in  an  otherwise  identical  radiation  field. 
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1.2.  Background 


1.2. 1.  Past  Methods  of  Experimental  Evaluation 

The  most  reliable  method  for  evaluating  this  topic  was  demonstrated  in  the 
APPLE  2  test  at  the  Nevada  Test  Site  during  OPERATION  TEAPOT  in  1955. 
Immediately  after  a  29  kT  detonation  atop  a  500-foot  tower,  approximately  1000  soldiers 
in  an  armored  task  force  called  Task  Force  RAZOR  drove  to  within  820  meters  of  the 
blast  site  as  the  mushroom  cloud  formed  overhead.  Upon  reaching  a  threshold  of  1  rad 
per  hour  inside  the  lead  tanks,  the  task  force  turned  and  maneuvered  to  assault  a 
designated  nearby  objective.  Radiation  monitors  tracked  exposures,  both  outside  and 
inside  the  vehicles,  as  well  as  to  all  participating  personnel  [2], 

While  the  APPLE  2  test  conclusively  proved  that  the  1955  fleet  of  vehicles  could 
operate  in  the  vicinity  of  a  nuclear  attack,  it  was  never  repeated.  Although  it  is  possible  to 
extract  some  information  about  the  radiation  protection  provided  by  the  vehicles  from  the 
test  data,  Task  Force  RAZOR  was  intended  to  demonstrate  capabilities  and  build 
confidence  of  a  tactical  formation,  not  to  specifically  evaluate  the  M48  tanks  and  M59 
annored  personnel  carriers  used.  After  the  Partial  Test  Ban  Treaty  of  1963  precluded 
further  tests  of  this  nature,  the  Army  developed  several  alternate  methods  for  evaluating 
RPFs,  which  it  used  for  the  subsequent  30  years.  These  approaches  included 
computational  modeling,  radioisotopes,  and  a  high-energy  X-ray  facility. 

As  the  U.S.  Army  fielded  new  vehicles  during  the  Cold  War,  each  chassis  was 
tested  and  evaluated  to  address  the  threat  of  tactical  nuclear  weapons  [3,  4,  15,  17,  18]. 
This  testing  evolved  with  improved  computing  power,  more  sophisticated  modeling 
approaches,  and  better  testing  apparatus.  The  evaluations  were  designed  to  provide  not 
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only  a  simple  scalar  RPF,  but  a  comprehensive  evaluation  that  would  advise  crews  and 
leadership  on  how  to  position  and  configure  their  vehicles  to  maximize  protection  against 
known  radiation  threats  [4,17,18]. 

Much  of  the  early  data  driving  gamma-based  RPF  assessments  used  radioisotopes 
with  actual  vehicles  and  radiation  detectors  [15].  These  were  often  measured  using  single 
isotopes,  most  commonly  60Co  [16].  The  complications  of  using  constantly  emitting, 
high-activity  sources  led  to  the  development  and  use  of  the  high-energy  X-ray  AURORA 
facility,  which  could  activate  or  shut  down  its  radiation  upon  command  [8].  Although  the 
still-active  Gamma  Irradiation  Facility  at  Sandia  National  Laboratories  was  designed  to 
accommodate  an  Ml  tank,  it  has  never  done  so. 

After  the  Cold  War  ended,  funding  and  interest  in  Anny  vehicle  RPFs  waned.  The 
AURORA  high-energy  photon  testing  facility,  which  evaluated  vehicle  gamma  protection 
factors  (GPFs),  ceased  operation  in  1994  [8].  An  effort  to  conduct  a  prompt  effects  test 
on  Army  vehicles  with  a  fast  burst  reactor  was  canceled  due  to  a  lack  of  funding.  There  is 
only  one  reactor  currently  available  which  can  irradiate  an  Anny  vehicle  with  a  prompt 
burst  of  fission  radiation,  the  Fast  Burst  Reactor  at  White  Sands  Missile  Range  [16].  The 
last  viable  alternative  for  a  prompt  fission  radiation  burst  was  the  Army  Pulsed  Radiation 
Facility,  which  was  shut  down  in  2003  and  demolished  in  2010  [1 1]. 

1.2.2.  Modeling  Efforts 

The  first  effort  to  model  RPFs  was  the  “Engineering  Method”  developed  by  the 
National  Bureau  of  Standards.  This  method  provided  semi-empirical  formulas  for 
transporting  gamma  radiation  through  common  materials  and  simple  geometries  [15]. 
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Although  this  approach  was  useful  for  designing  fallout  shelters,  it  lacked  the  flexibility 
needed  to  model  the  complicated  shapes  and  materials  of  a  vehicle. 

As  computing  technology  improved,  two  new  methods  became  viable:  discrete 
ordinates  transport  codes  and  Monte  Carlo  methods.  The  discrete  ordinates  transport 
technique  analytically  solves  the  radiation  transport  problem  in  tenns  of  the  average 
particle  with  numeric  integration  techniques,  whereas  Monte  Carlo  methods  use  random 
sampling  of  individual  particles  to  model  the  probabilistic  nature  of  neutral  particle 
interactions  [2].  Monte  Carlo  methods  emerged  as  the  better  modeling  option  in  the 
1980s.  A  crude  discrete  ordinates  transport  code  is  described  in  Chapter  3  and 
demonstrated  in  Chapter  5.  Monte  Carlo  methods  will  be  described  in  the  Chapter  2  and 
discussed  throughout  the  rest  of  this  paper. 

Although  these  efforts  have  been  able  to  model  parts  of  the  radiation  transport 
problem,  the  recent  release  of  MCNP6  presents  an  opportunity  to  model  all  aspects  of 
radiation  transport  in  a  single  program.  In  response  to  a  request  from  the  United  States 
Anny  Nuclear  and  Combating  Weapons  of  Mass  Destruction  Agency  (USANCA),  and  in 
coordination  with  the  Defense  Threat  Reduction  Agency  (DTRA),  a  new  research  effort 
began  in  2012.  This  effort’s  goal  is  to  produce  a  method  for  generating  reliable  radiation 
protection  factors  for  current  Anny  vehicles  [4,  14]. 

The  first  phase  of  this  project  was  completed  in  the  spring  of  2014,  when  the 
neutron  protection  factor  of  a  simple  geometric  surrogate,  modeled  using  the  Monte 
Carlo  Neutral  Particle  6  code  (MCNP6),  was  experimentally  validated.  The  simulation 
matched  experimental  measurements  with  a  confidence  level  of  96%  [4],  This 
investigation  continued  the  same  line  of  effort. 
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1.3.  Purpose 


The  goal  of  this  work  was  to  determine  the  suitability  for  MCNP6  as  a  model  for 
detennining  gamma  radiation  protection  factors  for  United  States  Anny  vehicles,  using  a 
simple  geometric  surrogate.  This  work  included  three  supporting  objectives.  The  first  was 
to  verify  the  suitability  of  MCNP6  by  using  other  modeling  approaches,  such  as  discrete 
ordinates  and  simplified  Monte  Carlo  methods.  This  verification  included  quantitative 
comparisons  and  qualitative  evaluations. 

The  second  objective  was  to  validate  MCNP6.  The  validation  used  a  cubic,  steel- 
walled  enclosure  with  a  22x22x22  inch  internal  cavity  as  a  surrogate  for  an  Anny 
vehicle.  The  validation  compared  experimental  data  from  a  detector  inside  and  outside 
the  enclosure  to  a  simulation  of  the  experiment  in  MCNP6.  The  threshold  for  successful 
validation  was  for  the  simulation  to  deliver  a  GPF  that  is  within  5%  of  the 
experimentally-derived  GPF  with  95%  confidence. 

The  final  objective  was  to  examine  the  influence  of  different  factors  on  the 
modeled  and  experimental  GPFs.  These  factors  included  using  different  photon  energies, 
different  geometric  configurations,  variance  reduction  techniques,  spectral  analysis  and 
processing  techniques,  and  different  model  fidelities.  This  objective  was  a  broad  effort 
that  included  physical  experimentation,  modeling,  and  data  processing  techniques  in 
order  to  analyze  these  influences  and  suggest  good  practices  for  future  research. 
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2.  Theory 


This  chapter  discusses  the  theory  of  gamma  radiation  transport  that  supports  this 
research,  methods  for  computing  dose  and  radiation  protection  factors,  the  theory  behind 
Monte  Carlo  methods,  information  about  MCNP6,  and  an  overview  of  alternate  models. 

2.1.  High-Energy  Photon  Interactions  with  Matter 

Individual  high-energy  photons  may  interact  in  several  ways  with  matter.  The 
location  and  type  of  interaction  are  driven  by  probabilities  known  as  cross  sections.  There 
are  two  categories  of  cross  sections:  microscopic  and  macroscopic. 

A  microscopic  cross  section  describes  the  probability  of  interaction  associated 
with  a  single  atom.  This  tenn  is  essential  for  determining  the  mean  free  path,  or  average 
distance  to  an  interaction  for  photons  of  a  given  energy.  A  macroscopic  cross  section 
sums  up  these  probabilities  for  all  types  of  atoms  in  a  material.  In  mathematical  tenns  [6], 

rx=i>i^)’ 

i=l 

where  £  is  the  macroscopic  cross  section,  crus  the  microscopic  cross  section  of  an 
atom,  p.  is  the  density  of  atoms  of  that  type  in  the  material,  and  n  is  the  number  of 
different  isotopes  in  the  material. 

The  microscopic  cross  section  can  be  further  broken  down  into  types  of 
interactions.  Gamma  rays  tend  to  interact  with  matter  through  three  well-understood 
mechanisms:  photoelectric  absorption,  pair  production,  and  Compton  scattering  [5,6]. 
The  total  microscopic  cross  section  is  a  simple  summation  of  the  microscopic  cross 
sections  of  each  type  of  interaction,  or  [6] 
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where  atot  is  the  total  microscopic  cross  section,  crpp  is  the  microscopic  cross  section  for 
pair  production,  a  is  the  microscopic  cross  section  for  photoelectric  absorption,  and 

pe 

acs  is  the  microscopic  cross  section  for  Compton  scattering.  The  microscopic  cross 

sections  of  these  interactions  depend  upon  the  number  of  protons  in  the  material  and  the 
energy  of  the  interacting  photon,  as  show  in  Figure  1 . 


Figure  1.  This  chart  shows  the  dominant  regions  of  each  type  of  interaction  based  upon  the  number 
of  protons  (Z)  in  the  material  and  the  energy  of  the  incident  photon.  The  lines  approximate  values  at 
which  neighboring  interaction  regions  are  equally  likely  [5]. 


In  photoelectric  absorption,  a  photon  interacts  with  an  absorbing  atom.  Some  of 
the  energy  of  the  photon  overcomes  the  binding  energy  of  one  of  the  atom’s  electrons, 
which  releases  that  electron.  The  rest  of  the  interacting  photon’s  energy  provides  the 
released  electron  with  kinetic  energy.  Since  charged  particles  have  mean  free  paths  far 
shorter  than  that  of  equivalent-energy  photons  [7],  secondary  electrons  from  this 
interaction  are  unlikely  to  travel  through  an  annor  plate.  Thus,  unless  the  photoelectric 
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absorption  happens  at  or  near  an  inner  surface  of  a  shield,  it  is  unlikely  to  transport  much 
energy  beyond  the  shield. 

In  a  pair  production  event,  a  high-energy  photon  interacts  with  an  atom’s  nucleus 
and  produces  two  beta  particles  of  opposite  charge.  In  order  to  balance  the  energy  of  this 
interaction,  the  incident  particle  must  have  energy  greater  than  or  equal  to  the  summed 
mass  energies  of  both  beta  particles.  Figure  1  shows  that,  for  the  heaviest  elements,  pair 
production  becomes  a  dominant  interaction  mechanism  for  photons  with  energies  over 
approximately  5  MeV.  Since  iron  has  26  protons  per  atom,  pair  production  becomes  a 
dominant  interaction  at  10  MeV,  which  is  far  more  energy  than  nearly  every  photon 
produced  in  a  nuclear  weapon  explosion  [2],  Since  the  Anny’s  interest  in  this  program  is 
driven  by  the  unique  effects  of  nuclear  weapons,  there  is  little  need  to  evaluate  photons 
with  energies  at  and  above  10  MeV. 

Compton  scattering  is  the  primary  interaction  of  interest  for  this  work.  In  a 
Compton  interaction,  a  photon  interacts  with  an  electron  and  deflects,  imparting  some  of 
its  energy  to  the  electron.  In  an  iron  barrier,  the  electron’s  mean  free  path  is  so  short  as  to 
make  its  effects  negligible  for  interactions  that  are  not  at  or  near  the  inner  surface  of  the 
barrier  [7],  but  the  photon  is  much  more  likely  to  continue  through  the  iron  barrier 
without  further  interactions.  The  amount  of  deflection  that  the  photon  experiences  is  a 
probabilistic  outcome  described  by  the  well-known  Klein-Nishina  distribution  [5,6], 

AP(Ef#f  p(E7.«)  +  P(E7,«f 1  -  1  +  cos2 (o) 


da 

cLC 
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(3) 


where  a  is  the  fine  structure  constant,  -  is  the  solid  angle  cross  section,  and  P  E  .0) 

d Q  V  7’  / 

is  the  ratio  of  input  and  output  photon  energies  given  by  the  Compton  scattering  equation 
[5,6], 


P 


(4) 


where  E,  is  the  energy  of  the  incident  photon  and  0  is  the  scattering  angle.  The  energy 

dependence  of  the  distribution  is  demonstrated  in  Figure  2  for  a  few  sampled  energies. 
Since  most  of  the  photons  of  interest  have  energies  between  100  keV  and  4  MeV,  the 
scattering  angles  will  tend  to  be  small,  so  that  the  scattered  photon  leaving  the  interaction 
will  have  most  of  the  energy  of  the  incident  photon. 


9 


90° 


Figure  2.  Polar  plot  of  Klein  Nishina  distributions  at  10  keV,  100  keV,  and  1  MeV  [13]. 


2.2.  Gamma  Protection  Factor 

In  its  simplest  fonn,  a  radiation  protection  factor  (RPF)  is  a  ratio  comparing  the 
effects  of  ionizing  radiation  in  a  protected  state  to  the  effects  of  ionizing  radiation  in  an 
unprotected  state  [14,  15], 

___  Free  Field  Dose 

Rrh  = - .  (5) 

Shielded  Dose 
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Since  charged-particle  radiation  at  energies  associated  with  nuclear  weapons  does  not 
have  any  significant  ability  to  penetrate  solid  matter  [5],  this  generalized  RPF  can  ignore 
the  influence  of  alpha  and  beta  radiation  as  external  sources.  However,  it  is  important  to 
note  that  photons  which  Compton  scatter,  or  interact  by  photoelectric  absorption,  at  or 
near  the  inner  surface  of  the  iron  barrier  may  release  energetic  electrons  that  subsequently 
contribute  to  the  dose  of  vehicle  occupants  as  a  secondary  source  of  ionizing  radiation.  If 
these  electrons  transport  significant  amounts  of  radiation  out  of  the  barrier,  then  they  will 
create  bremsstrahlung  radiation  and  increase  the  low-energy  photons  measured  in 
experimental  observations  above  what  is  predicted  by  a  model  that  does  not  include  these 
secondary  electrons.  Either  way,  these  electrons  are  a  secondary  product  of  gamma  or 
neutron  radiation,  so  the  RPF  formula  can  therefore  be  broken  down  into  a  sum  of  a 
neutron  and  a  gamma  protection  factor  [4,14,15], 


RPF  =  NPF  +  GPF, 


where 


NPF 


F  reeF  ieldNeutronDose 
ShieldedNeutronDose 


(6) 


(7) 


and 


GPF 


FreeFieldGammaDose 
ShieldedGammaDose  ’ 


(8) 


where  GPF  is  the  gamma  protection  factor  and  NPF  is  the  neutron  protection  factor.  How 
then,  should  we  define  the  free  field  and  shielded  gamma  doses?  In  general,  there  are 
three  approaches:  exposure,  dose,  and  dose  equivalent. 
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Exposure  is  defined  by  the  charge  created  within  a  mass  of  air  by  photon  radiation 


and  is  measured  in  Roentgens  (R),  which  are  defined  as  1R  =  2.58  *  10  4  0/  .  This 

/kg 

unit  of  exposure  is  often  used  in  the  study  of  radiation  protection  [5].  However,  this 
depends  on  cross  section  data  for  air,  and  therefore  must  be  adjusted  to  reflect  the 
differences  in  density  and  cross  section  between  air  and  a  target  material  of  interest,  such 
as  soft  tissue.  The  exposure  rate  can  be  computed  using  discrete  intervals  of  a  full 
spectrum,  using  [10] 
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where  X[E  is  the  photon-energy-dependent  exposure  rate  associated  with  a 

monoenergetic  flux,  (j)  is  the  flux  of  the  photons,  ph  is  a  number  of  photons,  Ephotong  is 
the  energy  of  the  photons,  p,  .  is  the  cross  section  associated  with  photons  of  that  energy 
in  air,  and  p  jr  is  the  density  of  air.  In  order  to  convert  the  exposure  rate  into  a  dose  rate, 
we  use  [10], 

D[Eph„o„J  =  X[EphotoJ  x  34  x  ,  (10) 

^ air  /  P air 

where  DI  E  ,  is  the  dose  rate  in  units  of  gray  per  second.  Conversions  must  be  done 

L  photons-1 

in  discrete  energies,  so  this  method  is  typically  applied  with  exposure  rate  constants  for 
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the  characteristic  energies  of  specific  radioisotopes  [5,9].  Given  that  both  prompt  and 
delayed  effects  of  nuclear  weapons  emit  broad  spectra  of  high-energy  photons  [2],  and 
given  that  MCNP6  does  not  have  a  tally  mechanism  for  exposure,  this  approach  can  be 
challenging  to  simulate,  especially  for  high-resolution  simulations  with  large  numbers  of 
bins. 


Another  method  is  to  evaluate  the  absorbed  dose  for  a  protected  body  against  the 
absorbed  dose  for  an  unprotected  body.  Dose  measures  the  total  energy  deposited  per  unit 
mass,  and  it  is  measured  in  units  of  gray  (Gy)  and  rad,  defined  as 

100  rad  =  1  Gy  =  1  ^  [2,5,6].  Dose  depends  upon  a  wide  variety  of  variables, 

including  the  geometry  of  the  object  of  study.  Since  human  bodies  inside  vehicles  may 
come  in  a  wide  variety  of  sizes  and  poses,  and  since  detectors  tend  to  have  different  cross 
sections  than  human  tissues,  a  sound  method  for  evaluating  GPFs  is  to  first  validate 
MCNP6  by  dose  depositions  in  experimental  and  simulated  detectors.  With  this 
validation  complete,  MCNP6  can  then  simulate  the  dose  deposited  in  a  volume  of 
surrogate  material  to  produce  a  reliable  human  body-based  GPF. 

Another  approach  is  to  compare  kenna  to  protected  and  unprotected  targets. 
Kenna  refers  to  the  total  kinetic  energy  deposited  by  ionizing  radiation  per  unit  mass 
[5,6,9].  Although  it  is  given  in  the  same  units  as  dose,  kenna  tends  to  be  higher  than  the 
actual  energy  imparted  to  a  material.  Some  of  the  energy  imparted  is  subsequently 
released  by  Bremsstrahlung  radiation,  which  can  escape  from  the  volume  of  interest 
[5,6,9].  This  feature  may  make  it  a  good  choice  for  a  worst-case  evaluation  of  a  body  or 
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object  in  MCNP6,  but  it  is  difficult  to  evaluate  experimentally,  since  all  secondary 
electrons  must  be  captured. 

2.3.  Other  Considerations 

Some  types  of  ionizing  radiation  have  more  adverse  health  effects  than  others  due  to 
high  coefficients  of  linear  energy  transfer.  This  is  particularly  true  of  heavy  charged 
particles  [9].  For  high-energy  photons,  however,  this  is  not  a  significant  factor,  since  the 
International  Commission  on  Radiation  Protection  (ICRP)  recommends  that  all  energies 
of  photons  be  assigned  a  weighting  factor  of  1,  the  lowest  weighting  factor  [5], 

Although  the  primary  penetrating  fonn  of  ionizing  radiation  for  this  study  is  the  high- 
energy  photon,  these  photons  deposit  their  energies  by  creating,  releasing,  and  exciting 
electrons.  Although  the  iron  shielding  is  likely  to  stop  most  of  these  electrons  due  to  their 
low  mean  free  path  in  iron,  any  electrons  released  or  excited  at  or  near  the  inner  surface 
of  the  iron  box  may  contribute  a  significant  amount  of  ionizing  flux  to  the  shielded 
environment.  The  ICRP  recommends  that  electron  radiation  be  assigned  a  weighting 
factor  of  1  [5]. 

Modem  Anny  vehicles  use  far  more  electronics  than  their  Cold  War  era  counterparts. 
This  introduces  the  question  of  how  vulnerable  the  vehicles’  systems  are  to  ionizing 
radiation.  This  may  be  a  valid  question  in  the  future,  but  most  modem  commercial 
electronics  require  50  Gy  or  higher  to  show  significant  degradation,  which  is 
approximately  25  lethal  doses  [10]. 
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2.4.  Monte  Carlo  Methods 


A  Monte  Carlo  approach  to  modeling  particle  transport  produces  an  approximation  of 
actual  data  by  tracking  individual  particles.  It  accomplishes  this  through  the  arrangement 
of  virtual  geometric  configurations  of  materials  to  produce  statistically  significant  results. 

This  approach  begins  by  randomly  sampling  an  energy,  direction,  and  distance  a 
particle  will  travel  until  it  interacts.  The  energy  depends  on  the  source  used.  The  direction 
is  a  random  sampling  of  an  even  distribution  of  angles.  The  distance  to  an  interaction  can 
be  sampled  by  the  equation 


dis  tan  ce 


—  ln(Rand(0, 1)) 

M2) 


(11) 


where  Rand(0,l)  is  a  random  real  number  between  0  and  1  [1].  This  produces  a 
probability  distribution  such  that  half  of  the  sampled  interaction  distance  is  between  one 
and  zero,  and  half  is  between  one  and  infinity. 

The  next  step  in  the  Monte  Carlo  process  is  to  transport  the  particle  in  its 
trajectory  until  it  reaches  its  interaction  distance.  This  depends  on  the  particle’s  mean  free 
path  length  in  whatever  medium  it  happens  to  pass.  The  mean  free  path  length  is 
detennined  by  the  particle’s  energy,  as  well  as  the  medium’s  density  and  experimentally 
derived  attenuation  data  [1], 

If  the  particle  crosses  into  a  new  material  before  its  interaction  distance  is 
exhausted,  a  new  interaction  distance  is  detennined.  Otherwise,  the  particle  interacts.  The 
type  and  outcome  of  the  interaction  are  energy-dependent;  see  section  2.1  for  more 
information.  Secondary  particles  are  generated  from  this  interaction,  and  these  particles 
are  tracked  by  this  same  methodology  until  they  meet  kill  criteria,  either  by  exiting  the 
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volume  of  interest  or  by  losing  so  much  energy  that  they  are  no  longer  a  significant 
source  of  ionization  [1], 

2.5.  Monte  Carlo  n-P article  Code  6  (MCNP6) 

Developed  by  Los  Alamos  National  Laboratories,  MCNP6  was  released  in  2013  [1], 
As  the  successor  to  a  series  of  programs  dating  back  to  the  first  Monte  Carlo  modeling 
efforts  in  the  1940s,  MCNP6  has  a  number  of  useful  features  that  go  far  beyond  neutral 
particles  to  include  36  types  of  particles  and  an  enonnous  library  of  interaction 
mechanisms  and  cross  sections  [1].  Its  predecessors,  MCNP5  and  MCNPX,  used  Monte 
Carlo  methods  to  simulate  different  types  of  particles,  and  were  developed  in  parallel  as 
separate  simulation  packages.  MCNP5  could  not  transport  charged  particles,  while 
MCNPX  could  not  transport  neutral  particles.  MCNP6  merged  the  capabilities  of  both 
programs  into  a  single  simulation  platform. 

User  inputted  text  files  drive  MCNP6  in  a  series  of  lines  known  as  ‘cards.’  The  first 
are  the  cell  cards,  which  describe  the  geometric  spaces  within  the  simulation,  their 
densities,  and  how  they  are  positioned.  The  second  are  the  surface  cards,  which  describe 
the  cell  boundaries.  Material  data  cards  define  the  proportions  of  elements  and  isotopes  in 
each  cell.  Other  data  cards  specify  the  particles  of  interest,  source  information,  types  of 
tallies,  bin  structure,  criticality  codes,  variance  reduction  techniques,  and  many  other 
options  [1].  Examples  of  MCNP6  input  files  are  in  Appendix  B. 
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2. 6.  Potential  Sources  of  Error 


Even  with  the  right  input  file,  MCNP6  does  not  account  for  dead  time  in  a  detector, 
the  relative  efficiency  of  a  detector,  or  the  response  function  associated  with  a  detector. 
Additionally,  if  the  input  deck  of  MCNP6  is  not  configured  to  track  electrons,  it  may 
ignore  a  significant  source  of  dose  within  a  protected  region. 

The  dead  time  in  a  detector  is  highly  dependent  upon  the  type  of  detector  and  the 
electronic  system  associated  with  that  detector.  The  Nyquist  criterion  requires  a  sampling 
rate  twice  that  of  the  measured  signal  in  order  to  resolve  the  signal.  If  the  system  cannot 
sample  fast  enough  to  establish  the  Nyquist  rate  for  two  or  more  photons,  then  the  system 
fails  to  record  the  extra  photons  [5].  This  is  a  significant  concern  when  using  the  same 
source  and  detector  in  both  a  shielded  and  unshielded  configuration.  An  unshielded 
configuration  with  a  strong  enough  source  could  lead  to  saturation  effects,  which  would 
degrade  the  utility  of  the  data  by  making  the  unshielded  activity  appear  lower  than  it  is. 
Timing  the  response  function  and  multiplying  it  by  the  number  of  counts  recorded  in  a 
given  operation  can  yield  the  total  dead  time  for  a  particular  run.  Subtracting  the  dead 
time  from  the  measurement  time  detennines  a  more  reliable  interaction  rate. 

Another  potential  source  of  error  is  electronic  system  noise.  If  the  signal-to-noise 
ratio  is  not  sufficiently  high,  then  data  is  lost  in  the  random  fluctuations  of  the  system. 
This  could  be  a  significant  problem  for  low-activity  scenarios  such  as  a  detector 
surrounded  by  the  iron  box  as  in  this  work. 

No  detector  is  perfect,  and  a  flux  of  monoenergetic  photons  into  a  spectroscopic 
instrument  will  register  as  a  sloped  peak,  not  a  single-bin  line.  This  detector  response  is 
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an  inherent  function  of  electronic  noise,  power  fluctuations  and  the  experimental 
tolerance  of  the  associated  components.  Since  MCNP6  tracks  actual  particle  behavior,  it 
will  put  all  particles  of  a  monoenergetic  flux  into  a  single  energy  bin.  Thus,  the  detector 
response  must  be  emulated  if  one  wishes  to  use  MCNP6  to  predict  actual  detector 
spectra.  The  full  width  and  half  maximum  of  a  signal  shows  the  difference  between  the 
independent  variable  at  half  the  maximum  height  of  the  dependent  variable  [5].  By 
measuring  the  full-width  at  half-max  of  a  known  peak,  it  is  possible  to  synthesize  a  point 
spread  function  which,  when  convoluted  with  a  MCNP6  generated  spectrum,  yields  a 
prediction  of  actual  experimental  data. 

Electrons  released  from  the  barrier  or  the  air  within  a  shielded  enclosure  may  pose 
another  issue.  Compton  interactions  are  the  primary  absorption  mechanism  of  interest  for 
transport  in  this  work,  but  all  three  types  of  high-energy  photon  interactions  can  release 
electrons  from  the  inner  surface  of  the  iron  box.  Although  these  electrons  are  not  as 
penetrating  as  the  down-scattered  photons,  they  may  contribute  to  dose. 


2. 7.  Hypothesis  Testing  and  Confidence  Intervals 

This  work  used  a  null  hypothesis  evaluation  to  detennine  if  the  GPF  produced  by 
MCNP6  was  or  was  not  valid.  The  initial  research  hypothesis  of  this  thesis  was  that  an 
MCNP6  model  of  a  shielded  and  unshielded  detector  subjected  to  the  same  external 
radiation  produces  a  GPF  within  5%  of  the  GPF  detennined  by  real  experiments  of  the 
same  configuration.  This  can  be  stated  in  mathematical  terms  as 


ModeledGPF  -  RealGPF 


<  .05  *  RealGPF. 


(12) 
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The  null  hypothesis  of  this  experiment  is  that  the  MCNP6  model  of  a  shielded  and 
unshielded  detector  subjected  to  the  same  external  radiation  produces  a  GPF  that  is  more 
than  5%  different  from  the  corresponding  GPF  determined  by  real  experiments  of  the 
same  configuration.  This  can  be  stated  in  mathematical  tenns  as 


ModeledGPF  -  RealGPF 


>  .05  *  RealGPF. 


(13) 


MCNP6  output  files  provide  the  variance  of  a  tallied  dose;  taking  the  square  root  of 
this  variance  yields  the  value  of  its  standard  distribution.  Assuming  that  the  GPF  is 
normally  distributed,  and  enough  points  are  sampled  by  both  the  experiment  and  the 
simulation,  a  two-tailed  test  with  a  p-value  of  0.05  should  disprove  either  the  hypothesis 
or  the  null  hypothesis  [11]. 
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3.  Experimental  Methods  and  Modeling 


This  chapter  describes  the  setup  and  execution  of  experiments  and  modeling 
simulations  to  verify  and  validate  MCNP6  as  a  model  for  producing  GPFs.  It  also 
explores  evaluation  techniques  and  alternate  experimental  configurations  for 
investigating  GPFs. 

3.1  Overview  of  the  Plan 

For  the  work’s  first  objective  of  verifying  MCNP6  as  a  suitable  choice  for 
modeling  gamma  transport,  four  alternate  models  were  compared  to  MCNP6  results.  The 
alternate  models  fell  into  two  categories:  non-Monte  Carlo  techniques  and  Monte  Carlo 
techniques.  In  an  effort  to  include  alternate  programming  languages,  these  alternatives 
were  written  on  different  programming  platforms.  All  told,  this  work  attempted  to  model 
photon  transport  in  Mathematica,  MATLAB,  GEANT4  (in  C++),  and  MCNP6,  which  is 
written  in  FORTRAN.  The  relative  merits  of  the  programming  languages  are  not 
addressed;  rather,  the  overall  processing  time,  portability,  versatility,  and  features 
provided  a  basis  for  evaluating  and  comparing  the  techniques. 

The  second  objective  of  validating  GPFs  produced  by  MCNP6  consisted  of  a 
simple  experiment  in  which  a  single-peak  source  irradiated  a  detector  in  a  simple  cubic 
steel  enclosure  and  in  a  free  field  environment.  The  experiment  was  run  twice,  once  with 
a  Nal  detector,  and  once  with  a  HPGe  detector. 

The  final  objective  of  evaluating  the  parameters  behind  GPF  validation  consisted 
of  several  series  of  simulations  and  experiments.  The  first  series  substituted  three 
different  sources  to  establish  the  range  of  energies  that  could  be  validated.  The  second 
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was  a  series  of  simulations  investigating  the  importance  of  geometric  fidelity,  material 
data  fidelity,  and  track  counts  to  achieving  validation.  The  third  series  investigated 
MCNP6’s  ability  to  produce  valid  GPFs  for  thicker  and  more-complicated  shielding 
geometries. 


3.2  Verifying  MCNP6  vs  Non-Monte  Carlo  Methods 

This  section  describes  the  methodologies  used  to  evaluate  non-Monte  Carlo 
photon  transport  techniques  against  MCNP6.  Two  models  were  constructed  and 
evaluated:  attenuation  and  quadrature. 


3.2.1  Verifying  MCNP6  vs  Attenuation 

Although  analytical  methods  could  not  model  the  precise  behavior  of  photon 
transport  above  the  quantum  scale,  it  was  possible  to  model  some  simpler  aspects  of 
photon  transport  with  a  single  equation:  the  exponential  attenuation  of  a  monoenergetic 
photon  beam.  Figure  3  is  a  representation  of  the  scenario  modeled.  This  model  made 
many  simplifying  assumptions,  such  as  directing  all  photons  on  a  perpendicular  path 
through  an  iron  barrier  of  finite  thickness  but  infinite  width  and  height. 

For  this  work,  the  model  was  constructed  in  Mathematica  to  provide  a  fast  and 
reliable  analytical  check  of  MCNP6  peak  photon  transport.  The  governing  equation  for 
monoenergetic  photon  attenuation  in  matter  is 

F  =  F()e-/tpx,  (14) 
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where  F  is  the  flux  of  a  stream  of  monoenergetic  photons  at  a  given  depth  in  a  barrier  x, 


Fq  is  the  original  flux  before  entering  the  barrier,  p  is  the  microscopic  cross  section  for  a 
photon  of  that  energy  in  that  material,  and  p  is  the  density  of  the  material  [6]. 


Figure  3.  Diagram  of  the  Mathematica  model,  with  uniform  flux  into  and  out  of  an  iron  barrier. 


The  cross  section  p  depends  on  both  the  materials  and  the  energy  of  the  photons, 
and  is  provided  in  material  data  tables  by  the  National  Institute  of  Standards  and 
Technology  (NIST)  [7].  Since  most  photon  energies  cannot  be  expected  to  have  the  exact 
reference  energies  in  these  tables,  the  cross  sections  must  be  interpolated  using  log-log 
interpolation. 

This  analytical  model  can  be  directly  compared  to  MCNP6  by  running  two 
simulations  with  a  monoenergetic  source  and  a  detector:  one  with  an  iron  barrier  between 
them,  and  one  without.  The  ratio  of  the  peak  energies  deposited  in  the  detector  should 
coincide  with  the  ratio  of  the  analytical  model’s  attenuated  peak  photons  to  the  original 
flux. 


0in  (analytical)  0Bare  (MCNP6) 
fatten  (analytical)  </>Shielded  (MCNP6)  ’ 


(15) 
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where  is  the  flux  going  into  the  shield  in  the  simple  attenuation  model,  </>attcn  is  the 
attenuated  flux  coming  out  of  the  other  side  of  the  shield  in  the  simple  attenuation  model, 
0Bare  is  the  flux  going  into  the  surface  of  the  unshielded  detector,  and  <Aghielded  is  the  flux 
going  into  the  surface  of  the  shielded  detector. 

3.2.2  Verifying  MCNP6  vs  a  Quadrature  Model 

Although  the  analytical  model  provides  reliable  and  verifiable  results  for 
monoenergetic  photon  penetration  and  attenuation,  it  does  not  provide  information  about 
scattered  photons.  The  Mathematica  model  from  the  previous  subsection  was  expanded  to 
demonstrate  this  behavior  for  two  successive  Compton  scatter  events.  This  approach  used 
numerical  integration  (or  quadrature)  techniques  by  breaking  up  all  probabilistic 
components  of  an  interaction  into  discrete  groups,  sub-groups,  and  sub-sub-groups  as  the 
model  progresses  through  one  and  two  subsequent  Compton  scatterings.  Although  crude 
by  modern  standards,  this  is  the  essence  of  how  discrete  ordinates  codes  work. 

All  interactions  are  sampled  in  discrete  increments  of  a  single  degree,  with 
nonnalized  weighting  depending  on  that  angle’s  probability  as  determined  by  the  Klein- 
Nishina  distribution  (See  Figure  2  in  the  previous  chapter). 

For  the  purposes  of  this  model,  all  photons  that  are  attenuated  in  the  previous 
section  are  assumed  to  Compton  scatter  at  a  distance  of  one  mean  free  path  into  the 
barrier.  Since  the  photons  going  into  this  first  Compton  scattering  are  all  nonnal  to  the 
plane  of  the  barrier,  photons  with  the  same  scattering  angle  all  have  identical  path  lengths 
through  the  rest  of  the  iron  barrier.  These  photons  must  then  be  attenuated  along  these 
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path  lengths  to  produce  a  spectrum  of  photons  that  scatter  forward  once  and  pass  through 
the  rest  of  the  iron  barrier. 

Further  Compton  interactions  follow  a  similar  logic  to  the  first;  all  attenuated 
photons  from  the  first  scattering  are  assumed  to  scatter  again  at  a  distance  of  one  mean 
free  path  further  into  the  iron.  However,  this  introduces  an  additional  complicating  factor, 
in  that  these  photons  are  no  longer  moving  nonnal  to  the  plane  of  the  barrier.  Thus  the 
scattering  angle  from  this  second  interaction  must  be  revolved  about  the  axis  of  the 
inbound  photon  to  this  second  interaction.  In  addition  to  further  dividing  the  photons  into 
180  scattering  angle  groups,  these  angles  must  be  divided  into  another  360  angles  of 
revolution. 

The  angles  of  the  inbound  photon,  the  incremental  scattering,  and  the  revolution 
of  that  scattering  about  the  axis  of  the  inbound  photon  produce  a  new  angle  with  respect 
to  the  nonnal  vector  of  the  barrier  plane.  This  is  done  through  the  Law  of  Cosines  for 
Spherical  Triangles  [12], 


Figure  4.  A  representation  of  a  Compton  interaction.  N  is  an  arbitrary  reference  vector,  which  is 
aligned  with  the  axis  of  the  iron  barrier  for  the  slab  geometry  of  the  Mathematica  and  MATLAB 
simulations. 
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Cos  (e)  =  Cos  (e ')  Cos  (0O )  +  Cos  (w)  sin  [e0 )  sin  [e ') , 


(16) 


where  9  is  the  angle  of  the  scattered  photon  with  respect  to  the  normal  vector  of  the 
barrier  plane,  6 1  is  the  input  angle  with  respect  to  the  nonnal  vector  of  the  barrier  plane, 
9q  is  the  scattering  angle  of  the  Compton  interaction,  and  u>  is  the  angle  of  rotation 

applied  to  the  scattered  angle  about  the  vector  of  the  inbound  photon.  Figure  4  is  a 
representation  of  the  angles  involved. 

The  resulting  model  should,  given  enough  interactions,  provide  a  rough 
approximation  for  the  expected  flux  of  a  photon  stream  through  an  iron  barrier.  Of 
particular  interest  is  the  number  of  subsequent  Compton  scatterings  that  this  approach  can 
practically  model. 

The  results  of  this  model  were  compared  to  a  MCNP6  model  using  both 
quantitative  and  qualitative  criteria.  The  quantitative  criteria  included  processing  time 
and  fde  sizes.  Although  this  model  produced  spectra,  they  were  compared  to  MCNP6 
spectra  in  tenns  of  the  spectra’s  shape  and  features,  not  the  counts  associated  with  those 
features.  The  qualitative  criteria  included  a  subjective  assessment  of  the  flexibility  of  this 
technique,  the  complexity  of  writing  an  input  file,  and  the  applicability  of  this  approach 
to  the  problem. 
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3.3  Verifying  MCNP6  vs  Other  Monte  Carlo  Programs 


This  section  describes  two  alternate  Monte  Carlo  programs  and  how  they  were 
configured  to  solve  the  same  problem  as  MCNP6.  The  first  program  was  a  custom 
MATLAB  script.  The  second  was  a  GEANT4  simulation. 

3.3.1  Slab  Geometry  Monte  Carlo  Model  in  MATLAB 

The  purpose  of  this  subsection  of  the  work  is  to  develop  a  Monte  Carlo  script  in 
MATLAB  that  models  photon  transport  through  a  barrier  of  finite  thickness  but  infinite 
height  and  width,  the  same  geometric  assumptions  as  the  previous  section’s  Mathematica 
model.  Unlike  that  model,  however,  each  photon  was  tracked  individually  by  randomly 
sampling  the  interaction  distributions  through  as  many  Compton  interactions  as  it  took  to 
meet  a  kill  criteria.  Photons  tracks  began  at  a  point  source  and  ended  either  by  passing 
out  of  the  volume  of  interest  or  losing  too  much  energy  to  be  of  further  interest. 
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Figure  5.  MATLAB  monte  carlo  model  slab  diagram. 

Figure  5  shows  a  representation  of  the  MATLAB  Monte  Carlo  script’s  geometry. 
The  slab  geometry  thicknesses  used  in  the  MATLAB  model  closely  matched  the  linear 
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dimensions  of  the  baseline  experimental  setup  along  a  ray  going  out  from  the  source 
through  the  center  of  the  Nal  detector  in  the  center  of  the  iron  box.  This  was  modeled  as  a 
plane  source,  a  slab  of  air,  a  slab  of  iron,  a  slab  of  air,  and  a  plane  representing  the 
surface  of  the  Nal.  The  result  was  an  energy  spectrum  of  photon  energies  crossing  into 
the  Nal,  which  could  then  be  multiplied  by  the  appropriate  geometric  and  dose 
conversions  to  approximate  the  spectrum  of  energies  deposited  for  the  MCNP6  model. 


3.3.2  GEANT4  Model 

This  section  describes  how  GEANT4  was  applied  with  an  identical  setup  to 
MCNP6.  GEANT4  is  an  open-source  Monte  Carlo  simulation  program  originally 
developed  by  CERN.  While  it  has  similar  simulation  capabilities  to  MCNP6,  GEANT4  is 
written  in  C++  and  requires  the  user  to  write  considerably  more  code  to  define  a 
simulation. 

The  models  were  set  up  to  simulate  the  planned  validation  experiment  described 
in  the  next  section.  Both  simulations  consisted  of  identical  iron  boxes  with  air  inside  and 
outside,  as  well  as  identical  Nal  detector  cylinders  placed  in  the  middle  of  the  boxes. 
Both  simulations  ran  662  keV  gamma  rays  from  isotropic  point  sources  at  the  same 
relative  positions  outside  the  steel  boxes  for  an  identical  number  of  tracks;  these  were 
meant  to  simulate  a  radioisotope  sample  of  137Cs.  The  energies  deposited  within  the  Nal 
detectors  were  reported  as  histograms  with  the  same  number  of  bins. 

Output  spectra  from  the  custom  MATLAB  program,  GEANT4,  and  MCNP6  were 
compared  by  qualitative  and  quantitative  criteria.  The  quantitative  criteria  included  a 
direct  comparison  of  photon  energy  spectra  deposited  in  the  Nal  detectors,  processing 
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time,  relative  peak  heights,  and  file  sizes.  The  qualitative  criteria  included  program 
flexibility,  input  file  complexity,  portability,  and  adaptability  to  different  aspects  of  a 
GPF  investigation. 

3.4  Validating  MCNP6  in  a  Simple  Experiment 

Although  establishing  validation  was  the  second  supporting  objective  of  this 
work,  it  was  the  most  important.  This  primary  experiment  compared  the  shielded  and 
unshielded  spectra  and  dose  rates  of  two  detectors  subjected  to  the  same  radioactive 
source  at  the  same  distance.  The  experiment  was  designed  to  be  as  simple  for  MCNP6  to 
model  as  possible.  Later  work  would  evaluate  more  complicated  scenarios  to  challenge 
MCNP6  and  explore  the  limits  of  evaluation. 

Figure  6.  Simple  diagram  of  a  shielded  observation,  not  to  scale. 

3.4.1  Experimental  Equipment  and  Setup 

The  test  enclosure  consisted  of  a  cubic  aluminum  frame  with  an  assortment  of  1 8 
steel  plates  which  could  be  bolted  onto  the  frame  to  fonn  a  shielding  barrier.  Each 
member  of  the  aluminum  frame  had  a  2  inch  by  2  inch  cross  section  and  outer  edges  22 
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inches  long.  Inside  the  enclosure,  a  light  aluminum  scissor  jack  held  the  detector  such 
that  the  detection  crystal  was  centered  in  the  middle  of  the  enclosed  space,  1 1  inches 
from  all  faces.  Figure  7  shows  the  frame  with  the  jack  and  stand  in  their  relative 
positions. 


Figure  7.  The  enclosure’s  components.  Top  left  shows  the  frame  with  the  stand  and  detector  inside. 
Top  right  shows  the  thicknesses  of  the  plates.  Bottom  shows  the  various  plates;  note  their  differences 
in  length  and  width. 

The  plates  themselves  were  made  in  sets  of  three  for  each  face  of  the  enclosure. 
Two  plates  of  each  set  had  thicknesses  of  x/i  inch,  while  the  third  had  a  thickness  of  % 
inch.  They  were  cut  to  fonn  a  cubic  shape  when  all  18  plates  were  bolted  onto  the 
aluminum  frame.  For  the  initial  validation  experiment,  only  the  inch  plates  were  used. 
In  this  thinner  assembly,  the  plates  overlapped  at  the  edges  and  comers,  as  shown  in 
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Figure  8.  The  source  was  therefore  placed  in  a  face-centered  position,  3  inches  away  from 
the  inner  surface  of  the  near  steel  plate. 


Figure  8.  Scale  model  of  the  assembled  enclosure  with  only  the  (4  inch  plates  on  (left)  and  all  three 
plates  on  (center),  and  the  source  placement  for  the  initial  validation  experiment  (right). 


The  Nal  detection  system  was  set  up  according  to  Figure  9.  The  Nal  detector  was 
a  Saint-Gobain  model  3M3/3  with  an  ORTEC  model  226  photomultiplier  tube  base.  The 
preamplifier  was  an  ORTEC  model  113.  The  EG&G  ORTEC  DART  was  a  multi-role 
component,  as  it  both  provided  the  high-voltage  power  supply  and  served  to  amplify  and 
digitize  the  signal.  Gammavision  allowed  the  desktop  computer  to  control  data  collection, 
display  it  in  real  time,  and  output  it  into  spectra  files  for  subsequent  data  processing. 


Computer  with 
Gammavision 


Figure  9.  Diagram  of  experimental  setup  for  Nal  detector  system. 
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The  Gammavision  settings  for  the  Dart  were  established  in  early  calibrations  and 
retained  throughout  the  rest  of  the  work.  The  high  voltage  power  supply  was  set  to  900  V. 
The  TTL  was  in  shutdown  mode.  The  ADC  gate  was  off,  its  conversion  gain  was  set  to 
8192  bins  to  allow  maximum  resolution,  and  the  lower  level  discriminator  was  set  to  88. 
The  energy  calibration  in  Gammavision  was  established  by  assigning  a  value  of  662  keV 
to  the  center  bin  of  the  observed  peak  and  applying  a  simple  linear  equation  going 
through  the  origin  of  the  spectrum. 

The  HPGe  detector  system  setup  was  similar  in  most  respects  to  the  Nal  detector 
system  setup.  The  most  significant  difference  was  the  inclusion  of  a  cryogenic  Sterling 
cooling  system.  This  consisted  of  an  ORTEC  Cryosecure  Compressor™  Power 
Controller,  an  ORTEC  X-Cooler™  II,  and  an  ORTEC  Pop  Top  Model  GEM-50 195-P- 
LP.  Figure  10  shows  the  setup  configuration  scheme  and  the  detector. 
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Figure  10.  HPGe  Detector  System  Diagram  and  picture  of  the  detector. 
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The  HPGe  detector  was  slightly  larger  than  the  opening  in  the  enclosure’s  top 
plate.  For  the  initial  validation,  the  top  panel  was  removed  to  avoid  the  necessity  of 
disassembling  and  reassembling  the  detector.  Since  the  detector  had  a  3 -day  cooldown 
and  wann-up  time,  this  presented  potential  delay  for  later  experimental  observations. 
However,  this  became  unnecessary,  as  the  HPGe  was  not  used  after  the  initial  validation 
experiment  for  reasons  discussed  in  Section  5.4. 

1  77 

Cs  was  selected  as  the  source  for  the  initial  experiment  due  to  its  single 
662  keV  photon  emission  peak.  Simulations  in  the  previous  sections  indicated  that  Cs 
would  penetrate  the  steel  in  significant  quantities  for  its  maximum  peak  and  in  its 
downscattered  regions.  The  source  was  placed  3  inches  outside  the  inner  surface  of  the 
steel  barrier,  or  2  14  inches  from  the  outer  surface.  After  taking  a  shielded  observation, 
the  detector  and  stand  were  removed  from  the  enclosure  and  placed  on  a  nearby  table. 
With  the  source  14  inches  away  and  at  the  same  height  as  the  detector  crystal’s  center, 
another  observation  was  taken. 

Rather  than  setting  a  fixed  observation  time,  the  end  of  an  observation  was  driven 
by  a  subjective  assessment  of  whether  or  not  the  key  features  of  a  spectrum  had  resolved 
with  an  acceptable  signal-to-noise  ratio.  The  spectra  produced  were  normalized  to  their 
observation  times  in  later  data  processing  in  MATLAB,  which  summed  the  weighted 
dose  contributions  of  each  bin  to  produce  a  single  dose  value.  Dividing  the  free  field  dose 
by  the  shielded  dose  yielded  the  experimental  GPF. 
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3.4.1  Model  Setup  and  Validation  Criterion 


The  input  deck  for  the  initial  validation  simulation  is  in  Appendix  B.  In  order  to 
keep  the  model  as  simple  as  possible,  the  corner  overlaps  were  not  included.  The  only 
physical  volumes  were  the  air  outside  the  enclosure,  the  steel  walls  of  the  enclosure 
(modeled  as  a  single,  solid  piece  with  no  corner  or  edge  overlaps)  and  a  Nal  cylinder 
suspended  in  the  center  of  the  enclosure.  The  source  was  modeled  as  a  point  source. 

Modeled  dose  was  computed  with  the  F8  tally,  which  counted  photon  tracks 
which  deposit  energy  within  the  detector  volume.  A  series  of  energy  bins  in  10  keV 
increments  allowed  the  program  to  produce  a  spectrum,  and  multiplying  these  bins  by 
their  medium  energy  and  summing  up  the  resulting  product  provided  a  modeled  dose.  As 
discussed  in  Section  1.3,  the  validation  criterion  was  arbitrarily  set  as  computing  the 
modeled  GPF  to  be  within  5%  of  the  measured  GPF  with  95%  confidence  or  greater.  The 
confidence  was  computed  using  the  statistical  data  provided  in  the  MCNP6  output  file. 

As  a  check  to  ensure  that  the  model  was  simulating  actual  behavior  properly,  the 
experimental  and  simulated  spectra  were  reproduced  using  the  simulation  and 
experimental  results.  This  allowed  side-by-side  comparisons  to  confirm  or  reject  the 
model’s  results.  These  results  are  discussed  in  Chapter  5  and  the  full  range  of  spectra  are 
presented  in  Appendix  A. 

3.5  Evaluating  Influences  on  Modeled  and  Experimental  GPFs 

The  final  objective  of  this  work  involved  several  series  of  simulations  and 
experimental  observations.  The  driving  purpose  was  to  investigate  the  limits  of  the 
validation  that  could  be  achieved  and  the  effects  of  different  factors  upon  modeled  and 
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experimental  GPFs.  The  first  series  studied  the  effects  of  different  photon  energies  upon 
GPFs  and  the  range  of  photons  that  upheld  the  validation  criterion.  The  second  examined 
how  different  types  and  levels  of  model  fidelity  affected  the  output  of  MCNP6.  Finally, 
several  series  of  simulations  and  experimental  observations  investigated  MCNP6’s  ability 
to  handle  thicker  shielding,  higher  GPFs,  and  more  complicated  geometries. 

3.5.1  Evaluating  GPFs  with  Different  Sources 

The  first  series  of  the  evaluation  objective  examined  the  question,  “How  does 
validation  hold  up  with  other  photon  energies?”  This  series  retained  the  same  model  and 
experimental  settings  as  the  initial  validation,  with  one  important  difference.  In  place  of 
the  Cs  source,  three  alternate  sources  irradiated  the  detector  in  its  shielded  and  un¬ 
shielded  configurations. 

The  first  alternate  source  was  60Co,  which  emits  two  peaks  at  the  same  activity. 
The  first  peak  is  1173  keV,  and  the  second  is  at  1332  keV.  Both  of  these  peaks  are 
substantially  higher  in  energy  than  Cs’s  peak  of  662  keV.  The  other  two  sources  were 
57Co,  which  emits  photons  at  122  and  138  keV,  and  109Cd,  which  emits  photons  at  88 
keV.  These  peaks,  along  with  their  Compton  edges,  provided  insight  into  how  low  and 
high  this  experimental  setup  could  confirm  validation. 

3.5.2  Evaluating  Model  Resolution 

The  next  series  in  the  evaluation  objective  was  conducted  with  simulations  only, 
and  it  was  intended  to  investigate  the  limits  of  MCNP6  and  determine  where  it  breaks 
down.  These  limits  were  evaluated  by  varying  three  parameters  of  the  input  deck: 
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geometric  fidelity,  material  data  fidelity,  and  the  number  of  tracks  evaluated.  Each 
parameter  featured  three  adjustments  to  a  specific  part  of  the  MCNP6  input  deck. 

Geometric  fidelity  is  a  subjective  assessment  of  how  closely  the  cells  and  their 
bounding  surfaces  match  the  physical  reality  they  are  meant  to  represent  in  the  model. 
This  parameter  had  three  levels  of  fidelity.  The  high  fidelity  model  was  the  initial 
validation  model.  The  medium  fidelity  model  eliminated  the  air  outside  the  enclosure  and 
removed  the  enclosure’s  wall  opposite  to  the  source.  The  low-fidelity  wall  removed  all 
walls  except  the  one  immediately  between  the  source  and  the  detector,  and  removed  all 
air,  redefining  the  rest  of  the  problem  as  a  vacuum. 

Material  data  fidelity  was  a  subjective  assessment  of  how  closely  the  model 
materials  matched  their  chemical  and  isotopic  proportions.  The  high  fidelity  model 
included  all  elements  in  a  given  material  making  up  more  than  1%  of  that  material,  as 
well  as  all  isotopes  in  that  element  making  up  more  than  1%  of  that  isotope.  The  medium 
fidelity  model  simplified  this  to  the  most  common  isotope  of  each  element.  The  low 
fidelity  model  simplified  this  still  further,  so  that  each  material  was  represented  by  the 
most-common  isotope  of  its  most  common  element,  except  for  Nal,  which  is  an  even 
compound  of  sodium  and  iodine. 

The  number  of  tracks  varies  the  number  of  primary  particles  created,  as  well  as  all 
of  each  primary  particle’s  secondary  particles.  This  varied  by  powers  of  10  from  103  to 
10  .  The  95%  confidence  interval  of  each  simulation  was  calculated  from  the  relative 
error.  By  identifying  the  lowest  number  of  track  counts  that  could  meet  validation 
standards,  a  minimum  threshold  could  be  established.  This  helped  the  overall  work  by 
eliminating  unnecessary  processing  time. 
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3.5.1  Evaluating  the  Full-Thickness  Assembly 


The  final  series  of  the  evaluation  objective  challenged  MCNP6  to  model  thicker 
shielding  and  more-complicated  geometries.  It  also  explored  how  these  geometries 
affected  GPFs.  The  enclosure  was  reassembled  using  all  18  steel  plates.  Figure  1 1  shows 
how  the  combined  thickness  of  these  plates  was  1  %  inches,  which  eliminated  the  comer 
overlap  of  the  thinner  assembly,  producing  a  cubic  shield  with  defined  comers  and  edges. 


Figure  11.  Scale  model  of  the  fully-assembled  box 

Since  higher-energy  photons  tend  to  penetrate  steel  more  often  and  were  already 
validated  for  the  14  inch  thickness  assembly,  the  137Cs  and  60Co  sources  were  used  for  this 
series  of  observations.  The  sources  were  evaluated  separately  in  three  locations.  Figure 
12  shows  the  placement  scheme  for  the  three  source  positions.  The  first  was  face- 
centered,  at  the  same  distance  from  the  center  of  the  detector  as  in  the  original 
experiment,  14  inches.  The  second  was  edge-centered,  at  the  same  height  as  in  the 
original  experiment,  but  placed  outside  the  edge,  so  that  the  distances  from  the  detector 
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center  to  the  source  along  both  horizontal  axes  of  the  box  were  14  inches.  The  third  was 
outside  the  corner,  such  that  the  distance  from  the  detector  to  the  source  along  all  three 
axes  of  the  box  was  14  inches. 


Figure  12.  Source  placement  diagram  for  the  1  54  inch  assembly.  FC  is  the  face-centered  position,  EC 
is  the  edge-centered  position,  and  CN  is  the  corner  position. 
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The  final  observation  used  the  two  sources  simultaneously.  The  Cs  was  placed 
at  the  edge-centered  position.  The  higher  energies  emitted  by  the  60Co  were  expected  to 
penetrate  more  steel ,  so  the  60Co  was  placed  at  the  comer.  Using  the  calibrated  activities 
of  these  two  sources,  they  were  modeled  in  a  single  MCNP6  simulation. 
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4.  Data  Processing 


This  chapter  describes  several  data  processing  techniques  used  in  the  work. 
MATLAB  served  as  the  post-collection  processing  platform  via  several  custom-made 
scripts.  These  scripts  computed  dose  depositions  within  the  detector,  convoluted 
simulated  spectra  to  provide  a  convenient  comparison  between  simulation  and 
experimental  observations,  and  produced  graphs  of  these  spectra.  Additionally,  the 
MCNP6  input  decks  employed  two  methods  of  variance  reduction,  and  Gammavision’s 
settings  were  adjusted  to  reduce  experimental  error  by  accounting  and  correcting  for  dead 
time. 

4.1  Computing  Dose 

The  dose  deposited  within  the  simulated  Nal  detector  and  the  experimental 
detector  were  drastically  different.  MCNP6  nonnalized  all  tallies  to  the  total  number  of 
primary  particles  run,  whereas  the  experimental  observation  is  a  simple  tally  for  each  bin. 
As  a  result,  the  graphs  of  experimental  spectra  are  labeled  as  “Number  of  Interactions”  or 
“Counts  per  Second,”  while  the  graphs  of  simulated  spectra  are  labeld  as  “Relative 
Activity.”  The  gamma  protection  factor  is  a  unitless  ratio  that  is  computed  using  only 
simulated  spectra  or  only  experimental  spectra.  Therefore,  no  conversion  is  necessary  for 
the  tallies. 

Computing  dose  with  MCNP6  was  a  straightforward  process.  The  F8  tally  in  the 
output  fde  provided  the  energy  of  each  bin  and  its  relative  activity.  Multiplying  each 
energy  by  its  activity  yielded  a  dose  associated  with  that  bin.  These  bin-specific  doses 
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were  then  summed  up  to  yield  a  net  dose  deposited  within  the  simulated  detector.  This 
process  is  represented  by 

n 

Dose  =  ^(EnergyjX  Activity;),  (17) 

i=l 

where  n  is  the  total  number  of  bins  in  the  F8  tally.  Dividing  the  dose  for  a  free  field 
simulation  by  the  dose  of  the  associated  shielded  simulation  yields  the  simulated  GPF. 

The  experimental  spectrum  posed  a  complication,  as  each  observation  included 
background  radiation.  To  correct  for  this,  three  multi-day  observations  were  taken  with 
no  sources  present  with  the  detector  outside  the  enclosure,  inside  the  14  inch  assembly, 
and  inside  the  1  %  inch  assembly.  The  output  file  from  Gammavision  included  the 
observation  time  with  the  full  array  of  the  tallies  recorded  in  each  energy  bin.  This 
allowed  these  background  spectra  to  be  normalized  by  dividing  the  tally  array  by  the 
observation  time 

i  Background 

iMormalized  Background  = - .  (18) 

Observation  Time 

The  experimental  spectra  were  then  corrected  using  these  nonnalized  background 
spectra.  First,  the  experimental  spectra  were  normalized  according  to  their  observation 
times  using  the  same  process  as  the  background  spectra.  Subtracting  the  appropriate 
nonnalized  background  spectrum  from  the  normalized  experimental  observation  yielded 
a  corrected,  nonnalized  experimental  activity. 

Using  the  same  equation  as  the  MCNP6  dose  calculation,  a  MATLAB  script 
produced  a  dose  per  unit  time.  The  same  process  was  used  (with  different  backgrounds) 
for  shielded  and  free  field  observations.  The  unit  of  time  introduced  by  the  background 
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correction  was  therefore  irrelevant  to  the  overall  GPF,  as  these  units  divided  out  to  yield  a 
unitless  ratio. 

4.2  Spectral  Analysis  and  Processing 

The  validation  criterion  of  this  work  was  a  comparison  of  simulated  GPFs  and 
their  confidence  intervals  to  a  corresponding  experimental  GPF.  This  criterion  was 
incomplete,  however,  without  assurances  that  MCNP6  was  modeling  the  interactions  and 
spectral  features  identified  by  the  Nal  detector  system.  Such  an  assurance  required  the 
application  of  signal  processing  techniques. 

When  a  photon  deposits  energy  into  the  actual  Nal  detector,  it  takes  a  sequence  of 
subsequent  physical  interactions  to  tally  that  deposition  in  Gammavision.  Many  of  these 
processes  are  inherently  imprecise,  such  as  the  process  in  which  the  photoelectrons 
created  in  the  detector  are  multiplied  through  a  series  of  electric  fields.  The  tendency  of 
the  detection  system  to  spread  out  a  signal  is  described  as  its  response  or  response 
function.  The  Nal  detector’s  response  to  a  discrete  peak  is  a  Gaussian  distribution  about 
that  peak,  while  MCNP6  will  show  the  discrete  peak  with  no  Gaussian  distribution. 
Figure  13  compares  an  MCNP6  spectrum  to  an  experimental  spectrum  for  similar 
configurations  of  sources  and  detectors. 


40 


MCNP6  Energy  Spectrum 


Initial  Validation.  Shielded  Experiment 
Cs-137  Source  15  Inches  from  Nal  Detector 


Figure  13.  An  unprocessed  MCNP6  spectrum  and  an  unprocessed  experimental  spectrum.  The 
maximum  peak  of  the  MCNP6  spectrum  is  well  off  the  scale. 


Manipulating  the  simulation  to  emulate  the  response  of  the  detector  requires 

i  'in 

convolution,  which  is  described  in  Section  2.6.  Since  Cs  provides  a  single  662  keV 
peak,  its  maximum  peak  serves  as  a  convenient  source  for  evaluating  the  detector 
response.  The  full  width  at  half  maximum  was  measured.  This  value  was  calculated  in 
terms  of  the  total  number  of  counts  detected  within  the  peak  as 

FWHM  =  0.15  *  N,  (19) 


where  N  is  the  total  number  of  counts  within  the  peak.  Since  this  peak  follows  the  pattern 
of  a  Gaussian  distribution,  it  follows  that  the  standard  deviation  can  be  computed  as  [22] 


a 


FWHM 
2^2  ln(2)  ’ 


(20) 


where  FWHM  is  the  full  width  at  half  maximum  and  cris  the  standard  deviation.  To 
apply  this  to  a  discrete  digital  signal,  the  axis  of  the  signal  must  be  defined  as  a  linear 
space  vector  of  N  points  between  -N/2+1  and  N/2.  This  vector  is  defined  in  MATLAB  as 
X  using  the  linspace  function: 

X  =  linspace(-N/2+l,N/2,N).  (21) 
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The  point  spread  function  (PSF)  of  a  Gaussian  distribution  is  then  defined  as  [22] 


x2 

PSF  =  — =!^=e  2^2. 

\27UJ2 


(22) 


By  taking  the  Fourier  transform  of  this  point  spread  function,  the  Fourier  transform  of  the 
signal,  multiplying  them  together,  and  taking  the  inverse  Fourier  transform  of  the  result, 
MATLAB  generated  a  padded  array  containing  a  convoluted  spectrum  with  the  PSF 
applied  to  every  bin.  Once  this  array  was  trimmed  of  the  extra  0  value  elements  and 
calibrated  to  the  known  peak,  it  produced  a  simulated  spectrum  that  closely  matched  the 
experimental  observation,  as  shown  in  Figure  14. 


Initial  Validation,  Shielded  MCNP6  Simulation,  10  Million  Tracks 
,  Cs-137  Source  15  Inches  from  Nal  Detector 


Energy  of  Interaction  (keV) 


Figure  14.  Comparison  of  a  convoluted  simulated  spectrum  (left)  with  its  experimental  counterpart 
(right). 


4.3  Variance  Reduction 

Although  it  is  not  a  major  component  of  this  work’s  investigation,  variance 
reduction  is  an  important  aspect  of  MCNP6  modeling,  and  several  of  its  key  principles 
were  applied  to  the  simulations.  Variance  reduction  is  a  collection  of  processes  by  which 
the  statistical  resolution  of  data  collected  can  be  improved  while  reducing  the  overall 
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computation  time  of  the  simulation.  There  are  three  principle  methods  of  variance 
reduction  in  MCNP6. 

The  first  technique  of  variance  reduction  is  known  as  truncation.  Truncation  is  the 
practice  of  killing  particles  that  leave  the  geometry  of  interest.  In  this  work,  the  principle 
of  truncation  was  applied  by  reducing  the  surface  bounding  the  problem  to  be  a  meter 
away  from  the  enclosure’s  inner  edge  for  the  %  inch  assembly,  and  two  meters  away 
from  the  enclosure’s  inner  edge  for  the  1  %  inch  assembly. 

Russian  roulette  and  particle  splitting  involve  reducing  the  computation  time 
associated  with  regions  of  low  statistical  significance.  When  a  particle  enters  such  a 
region,  a  weighted  random  sampling  decides  whether  or  not  to  kill  the  particle  or 
continue  to  track  it.  This  is  known  as  Russian  roulette.  If  the  tracked  particle  passes  back 
into  a  region  of  interest,  it  is  split  into  several  particles,  each  of  which  are  then  tracked. 
This  is  known  as  particle  splitting.  In  MCNP6,  these  practices  are  done  by  adjusting  the 
importance  of  a  cell.  For  example,  in  the  1  Vi  inch  assembly  model,  a  concrete  room  was 
modeled  around  the  enclosure.  This  increased  the  simulation  run  time  on  the  laptop 
computer  used  from  88  seconds  to  547  seconds.  By  assigning  an  importance  of  .1  to  the 
concrete  cell,  and  an  importance  of  1  to  everything  else,  the  runtime  was  reduced  to  152 
seconds,  while  the  statistical  variance  in  the  detection  cell  was  unchanged. 

4.4  Eliminating  Dead  Time 

As  discussed  in  Chapter  2,  the  dead  time  of  a  detector  is  an  important  feature  in 
the  validation  that  is  not  modeled  in  MCNP6.  This  is  a  particularly  concerning  factor  for 
experiments  in  which  detectors  experience  very  high  dose  rates.  If  not  considered  and 
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corrected  for,  this  could  be  a  problem  for  unshielded  detectors.  Fortunately,  most  spectral 
software  keeps  track  of  the  dead  time,  allowing  the  flux  to  be  corrected  with  a  simple 
operation  for  each  energy  bin: 


Counts^  ,  , 

Corrected 


Counts^ 

Haw 


* _ TimeQbserved _ 

TimeObserved-TimeDead’ 


(23) 


Gammavision  was  set  to  apply  this  correction  to  all  experimental  spectra  taken. 
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5.  Results  and  Analysis 


This  chapter  presents  results  of  models,  experiments,  and  data  processing  to 
support  the  three  objectives  of  this  work:  verifying  MCNP6  as  a  good  choice  for  GPF 
computation  against  other  approaches  and  models,  validating  MCNP6  for  GPF 
computation  with  a  simple  experiment,  and  evaluating  the  influences  of  different  factors 
upon  validation  and  GPFs.  Results  supporting  the  objective  of  verifying  MCNP6  as  a 
valid  modeling  choice  are  presented  in  Sections  5.1  and  5.2.  Section  5.3  validates 
MCNP6  against  a  simple  experiment.  Sections  5.4  through  5.8  evaluate  the  influences  of 
several  factors  upon  validation  and  GPFs. 

5.1  Comparison  of  Analytical  and  Quadrature  Models  to  MCNP6 

The  Mathematica  model  was  set  to  have  an  initial  flux  of  100  keV  photons 
through  a  1  cm  iron  barrier.  The  analytical  portion  of  this  model  provided  a  means  for 
direct  quantitative  comparison  with  MCNP6.  For  the  Mathematica  model,  a  “Peak  GPF” 
was  computed  by  dividing  the  attenuated  flux  by  the  original  flux. 

In  MCNP6,  a  1  cm-thick  slab  was  placed  between  a  source  and  a  Nal  cylinder. 
Two  simulations  were  run,  one  with  the  barrier  and  one  without.  Both  simulations  tallied 
the  full-peak  photon  energy  depositions  inside  the  simulated  detector.  The  unshielded 
peak  activity  was  then  divided  by  the  shielded  peak  activity.  The  results  of  the  MCNP6 
simulation  and  the  analytical  simulation  are  in  Table  1.  The  95%  confidence  interval  for 
the  MCNP6  result  is  two  standard  deviations,  while  the  analytical  result  is  not  derived 
from  statistics  and  therefore  has  no  confidence  interval. 
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Table  1.  Peak  gamma  protection  factors  for  MCNP6  and  the  analytical  model. 


Peak 

GPF 

2  a 

Analytical 

1.783 

N/A 

MCNP6 

1.746 

0.0028 

%Difference 

2.12 

This  comparison  verified  MCNP6’s  peak  attenuation.  However,  it  provided  no 
information  on  the  scattered  photons.  The  scattered  photons  were  modeled  in 
Mathematica  through  two  successive  Compton  scatterings. 

All  of  the  photons  that  attenuated  through  the  barrier  were  assumed  to  Compton 
scatter.  The  first  scattering  was  modeled  as  happening  in  the  center  of  the  barrier,  and  the 
second  was  modeled  as  happening  3A  of  the  way  through  the  barrier.  Figure  15  shows  the 
binned  photons  escaping  after  the  first  two  scatterings.  The  second-scatter  escapes  are 
more  plentiful  than  the  first  due  to  the  fact  that  they  had  only  half  the  thickness  of  iron 
left  to  penetrate. 


Figure  15.  Binned  photon  energy  groups  that  pass  through  the  iron  barrier  after  one  scattering 

event  (red)  and  two  scattering  events  (blue). 
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The  spectrum  produced  by  the  first  two  scatterings  shows  some  interesting 
differences  between  the  first  and  second  scatterings.  In  particular,  the  second  scattering 
activity  peaks  at  a  lower  energy,  while  the  first  scattering  peak  is  essentially  at  the 
original  peak  energy.  This  demonstrates  the  beginnings  of  what  would  likely  develop  into 
an  escape  peak  from  the  iron.  Further  scatterings,  however,  proved  to  be  impractical  to 
model. 

The  first  scattering  required  Mathematica  to  calculate  the  probability  of  180 
degrees  of  scattering  angles  in  180  separate  groups.  These  groups  were  then  attenuated 
through  the  remaining  half  of  the  barrier.  For  the  second  scattering,  the  attenuated 
portions  of  these  180  groups  had  to  be  divided  again  with  180  subgroups  for  each  of  the 
original  180  groups.  Additionally,  the  scattering  angle  had  to  be  revolved  around  the  axis 
of  the  incident  photon  in  accordance  with  the  Law  of  Cosines  for  Spherical  Triangles 
described  in  Chapter  2.  This  required  each  subgroup  to  be  divided  into  360  sub¬ 
subgroups.  Each  sub-subgroup  then  had  to  be  attenuated  through  the  remaining  quarter  of 
the  barrier.  Any  further  scatterings  would  have  required  64,800  times  as  many 
calculations  as  the  one  before  it. 

This  geometric  rise  in  required  computing  became  impractical  beyond  two 
successive  scatterings.  Using  an  Intel  Core  i5-2520M  CPU  at  2.5GHz  with  4  GB  of  RAM 
and  no  other  active  programs,  Mathematica  took  150  minutes  to  generate  the  outcome  of 
two  successive  scatterings.  Another  scattering  would  increase  this  time  to  18.5  years. 
Although  this  program  could  be  improved  by  eliminating  groups  with  activities  below  a 
certain  threshold  and  using  fewer  groups  per  scattering,  the  geometric  expansion  of 
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computing  time  for  successive  Compton  scatters  would  remain.  Without  more  Compton 
scatters,  the  program  could  not  produce  a  comparable  spectrum  to  MCNP6.  Without 
similar  spectral  behavior  between  the  two  models,  there  was  no  basis  to  perform  a 
meaningful  quantitative  GPF  comparison.  These  results  were  useful  for  understanding 
some  of  the  phenomena  of  photon  interactions,  but  they  were  not  useful  for  modeling 
GPFs. 

Furthermore,  the  model  was  inflexible  and  required  extensive  retooling  to  change 
between  different  materials.  It  could  not  model  complex  geometries,  and  its  file  size  was 
1.1  megabytes,  as  compared  to  an  MCNP  input  file  of  15  kilobytes.  Finally,  Mathematica 
itself  is  expensive,  costing  $2,250  for  a  standard  government  license  [12].  By  contrast, 
MCNP6  is  available  upon  request  from  LANL  and  does  not  require  additional  taxpayer 
funding. 

5.2  Comparison  of  Other  Monte  Carlo  Models  to  MCNP 6 

Since  the  MATLAB  slab-geometry  Monte  Carlo  program  was  not  limited  to  a 
certain  number  of  successive  interactions,  it  produced  a  spectrum  that  was  more 
comparable  to  the  MCNP6  model  than  the  Mathematica  approach.  The  spectrum 
produced  was  a  binned  tally  of  the  energies  of  photons  crossing  a  2-dimensional  barrier 
between  two  iron  barriers.  Figure  16  compares  the  MATLAB  Monte  Carlo  simulation’s 
spectrum  to  an  MCNP6  tally  of  photons  crossing  the  surface  of  a  detector  within  an 
analogous  steel  enclosure. 
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MATLAB  Monte  Carlo  Simulation  Simulated  Shielded  Flux  into  HPGe 


Photon  Energy  (McV)  Energy  of  Photons  Crossing  Detector's  Surface  (keV) 


Figure  16.  MATLAB  Monte  Carlo  script  output,  simulating  the  photon  flux  between  two  iron 
barriers  (left)  and  MCNP6  simulation  of  photon  flux  into  a  detector  (right). 

Although  both  spectra  show  an  escape  peak  just  under  200  keV,  the  MATLAB 
model  also  shows  a  larger  peak  at  300  keV  and  a  high-energy  activity  that  rises  gradually 
to  the  peak.  Neither  of  these  behaviors  is  supported  by  theory.  The  cause  of  this  behavior 
was  not  determined,  but  it  was  likely  due  to  a  bug  in  the  transport  code  or  a  side-effect  of 
using  slab  geometry.  Without  similar  spectral  behavior  between  the  two  models,  there 
was  no  basis  to  perfonn  a  meaningful  quantitative  GPF  comparison. 

Qualitatively,  MATLAB  offered  some  significant  benefits,  including  easy 
debugging  features,  excellent  native  graphing  tools,  and  a  code  structure  that  was  short 
and  easy  to  examine  in  detail.  Further  work  might  well  be  able  to  improve  the  model. 
However,  even  with  its  current  issues  resolved,  the  slab  geometry  model  could  not 
account  for  complex  3 -dimensional  shapes  or  model  any  other  interactions  besides 
Compton  scattering.  The  model  would  also  require  significant  reprogramming  for 
different  layer  configurations,  materials,  and  tallies.  The  MATLAB  model  was  also  less 
efficient  than  MCNP6  for  large  track  numbers.  While  using  the  same  processor  with  no 
other  programs  running,  a  simulation  of  10  million  tracks  took  MCNP6  two  minutes, 
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while  the  MATLAB  code  required  two  hours.  Finally,  MATLAB  costs  $2,150  per 
government  or  commercial  license,  while  MCNP6  is  free  [21]. 

Of  the  modeling  options  examined,  GEANT4  comes  closest  to  MCNP6  in  overall 
performance.  Both  are  available  at  no  cost,  and  both  can  model  a  wide  variety  of  particles 
through  three  dimensional  geometries  and  complex  materials.  In  some  ways,  GEANT4  is 
superior  to  MCNP6.  For  example,  it  retains  individual  particle  track  information,  and  its 
native  visualization  allows  users  to  ensure  a  correct  geometric  configuration  and  view 
particle  tracks  directly,  as  seen  in  Figure  17. 


Figure  17.  GEANT4  visualization  of  the  model  with  some  demonstration  photon  tracks.  The  source  is 
outside  the  box,  on  a  center  line  with  the  nearest  face.  Inside  is  a  Nal  cylinder. 

The  identical  inputs  to  GEANT4  and  MCNP6  yielded  energy  deposition  spectra 
in  the  Nal  detector  volumes.  Both  spectra  were  plotted  in  67  energy  bins  10  keV  wide  as 
shown  in  Figure  18. 
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Figure  18.  Photon  energy  deposition  spectra  for  identical  inputs  to  MCNP6  (left)  and  GEANT4 
(right).  The  peak  spectra  at  662  keV  are  cut  off  in  order  to  show  down-scattered  photon  energy 
spectra. 


A  cursory  examination  of  these  spectra  highlights  their  similar  features.  Both 
models  produced  peaks  at  662  keV,  Compton  edges  at  460  keV,  and  build  up  regions  at 
the  same  photon  energies  around  200  keV.  The  overall  shapes  of  the  spectra  match  up 
well.  Figure  19  plots  the  difference  in  each  bin  from  GEANT5  to  MCNP6.  Only  one  bin 
differs  by  more  than  50%  from  its  counterpart  in  the  other  spectrum. 
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Energy  Spectrum  Difference  (GEANT4-MCNP6) 


Figure  19.  A  difference  plot  between  the  GEANT4  energy  deposition  spectrum  and  the  MCNP6 
energy  deposition  spectrum.  Ideally,  this  would  be  a  flat  line  at  zero. 

The  primary  discrepancy  between  these  spectra  was  the  peak  monoenergetic 
photons  that  passed  through  the  iron  without  interacting  and  deposited  all  their  energy 
inside  the  detector  volume.  GEANT4  showed  9,186  interactions,  while  MCNP6  showed 
only  5,589.  Additionally,  there  is  a  sudden  drop  off  in  the  GEANT4  spectrum’s  lower- 
energy  build-up  region  that  is  not  present  in  the  same  region  of  the  MCNP6  region. 

One  possible  discrepancy  between  the  two  codes  is  that  GEANT4  is  configured  to 
track  electrons,  while  MCNP6,  for  this  run,  was  not.  However,  this  explains  neither  the 
difference  in  peak  photon  counts  nor  the  different  build-up  region  shapes,  since  tracking 
secondary  electrons  would  result  in  GEANT4  having  more  activity  at  lower  energy 
levels,  and  not  less.  Another  possible  explanation  is  an  error  in  one  of  the  simulation  files 
for  GEANT4,  which  requires  significantly  more  programming  than  MCNP6. 
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It  is  difficult  to  assess  which  model  is  more  correct  by  inspection  of  the  spectra, 
since  no  detector  will  give  a  perfect,  100%  efficient  energy  collection  with  perfect 
resolution.  However,  it  is  possible  to  predict  the  actual  detector  spectrum  using  signal 
processing  techniques.  Using  the  detector’s  measured  full  width  half  maximum  of  a 
known  photopeak  for  a  radioactive  isotope,  it  is  possible  to  construct  a  system  response 
as  a  point  spread  function  for  the  detection  system.  Convoluting  the  simulated  spectrum 
with  this  point  spread  function  yields  a  close  approximation  of  a  Nal  detector  spectrum. 


MCNP6  Energy  Spectrum  (Convoluted) 


Energy  Deposited  (keV) 


GEANT4  Energy  Spectrum  (Convoluted) 


Figure  20.  A  comparison  of  energy  deposition  spectra  convoluted  with  a  point  spread  function  to 


simulate  an  actual  Nal  detector. 


Figure  20  shows  the  two  simulators’  convoluted  spectra.  Comparing  later 
experimental  spectra  to  these  convoluted  simulation  spectra  indicated  that  the  MCNP6 
simulation  was  a  better  representation  of  reality.  However,  for  the  purposes  of  this  work, 
either  model  might  produce  valid  GPFs  regardless  of  the  outcome  of  such  a  comparison, 
as  long  as  the  modeled  ratio  of  shielded  dose  to  unshielded  dose  is  consistent  with 
experimental  results. 
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5.3  Initial  Validation  Experiment:  Nal 


Comparing  dose  depositions  with  MCNP6  to  the  dose  detected  in  experimentally 
shielded  and  unshielded  configurations  yielded  GPFs  with  a  difference  of  3.4%  for  the 
Nal  detector.  The  high  number  of  MCNP6  tracks  and  large  number  of  interactions 
measured  reduced  the  variance  to  negligible  levels;  this  is  explained  in  detail  in  Section 
5.5.3.  Therefore,  MCNP6  met  the  validation  criteria  for  this  experiment  established  in 
Section  2.7. 

Table  1  presents  the  experimental  and  modeled  GPFs,  as  well  as  their  95% 
confidence  intervals,  defined  as  the  region  within  two  standard  deviations  of  the  GPF 
presented.  Since  GPFs  are  dose  ratios,  no  units  are  necessary  for  this  comparison. 

Table  2.  Initial  validation  results 


MCNP6 

Experimental 

GPF 

1.0225 

1.0587 

±2  a 

0.0016 

0.0053 

In  this  case,  a  simple  ratio  comparison  raises  some  logical  concerns  about  the 
reliability  of  this  evaluation.  After  all,  dividing  the  dose  deposited  from  one  free-field 
measurement  by  an  identical  measurement  would  yield  a  GPF  of  1,  which  is  almost 
within  the  validation  criterion.  Therefore,  while  the  GPF  comparison  met  the  validation 
threshold,  it  requires  some  supporting  spectral  comparisons. 

Figure  21  and  Figure  22  show  a  direct  comparison  between  spectra  created  from 
MCNP6  data  and  spectra  taken  experimentally.  The  consistency  of  their  peaks  and 
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Compton  edges  supports  the  GPF  validation.  However,  there  are  some  interesting,  if 
minor,  discrepancies. 


Initial  Validation,  Free  Field  MCNP  Simulation.  10  Million  Tracks 
Cs-137  Source  15  Inches  from  Nal  Detector 


Initial  Validation,  Free  Field  Experiment 
Cs-1 37  Source  1 5  Inches  from  Nal  Detector 


Energy  of  Interactions  (keV) 


SCO  TOO 


Figure  21.  Simulated  (left)  and  experimental  (right)  spectra  for  a  free-field  configuration.  The 
simulated  spectrum  has  been  convoluted  with  a  point-spread  function  to  generate  an  ideal  Nal 
response.  The  experimental  results  do  not  include  the  lowest  energy  bins,  as  these  tend  to  become 
saturated  with  electronic  noise. 


Although  these  results  met  expectations,  several  differences  bear  mentioning.  The 
free-field  spectra  of  the  MCNP6  and  experimental  results  show  identical  energy  values 
for  the  peak  and  Compton  edge,  but  the  experimental  results  show  some  build-up  at 
lower  energies,  as  well  as  a  very  low-energy  peak.  A  combination  of  background 
radiation  and  photon  interactions  with  the  nearby  concrete  and  low-energy 
phosphorescence  in  the  detector  are  the  likely  causes  of  these  discrepancies.  Since  these 
discrepancies  occur  in  the  low-energy  bins,  they  have  little  impact  on  the  overall  dose  or 
GPF. 

These  discrepancies  were  notably  absent  in  the  experimental  shielded  spectrum, 
which  closely  matched  the  simulated  spectrum’s  features  as  shown  in  Figure  22.  The  only 
notable  difference  was  a  lower-energy  build-up  region  with  noticeably  more  interactions 
for  the  experimental  results  than  the  simulated  results.  Since  this  discrepancy  paralleled 
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the  increased  low-energy  interactions  of  the  free  field  experimental  spectrum,  the  overall 
GPF  remained  consistent. 

Initial  Validation,  Shielded  MCNP6  Simulation,  10  Million  Tracks  Initial  Validation,  Shielded  Experiment 


Cs-137  Source  15  Inches  from  Nal  Detector  Cs-137  Source  15  Inches  from  Nal  Detector 


Energy  of  Interaction  (keV)  Energy  of  Interactions  (keV) 

Figure  22.  Simulated  (left)  and  experimental  (right)  spectra  for  a  shielded  configuration.  The 
simulated  spectrum  has  been  convoluted  with  a  point-spread  function  to  generate  an  ideal  Nal 
response.  The  experimental  results  do  not  include  the  lowest  energy  bins,  as  these  tend  to  become 
saturated  with  electronic  noise. 

The  limitations  of  the  equipment  itself  create  additional  sources  of  error.  The 
response  function  of  the  detector  spread  the  662  keV  single  energy  peak  more  than  the 
experimentally-derived  point  spread  function  used  to  convolute  the  model  results.  This  is 
likely  due  to  fluctuations  in  the  power  supply.  Additionally,  the  lowest  energy  bins  had  to 
be  truncated,  as  these  bins  tended  to  become  saturated  with  electronic  noise  and 
phosphorescence  from  the  Nal.  Additionally,  these  spectra  were  calibrated  only  with  the 
662  keV  energy  peak  with  the  assumption  that  the  Nal  detector  had  a  linear  response 
function.  A  more  thorough  multi-peak  calibration  could  match  the  energy  bins  somewhat 
more  precisely.  Although  both  of  these  could  have  been  corrected  with  additional 
equipment  and  a  more  extensive  series  of  calibration  runs,  neither  factor  prevented 
validation,  so  their  role  was  considered  negligible  for  the  purposes  of  this  study. 
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5.4  Initial  Validation  Experiment:  HPGe 


The  energy-dependent  efficiency  of  the  gennanium  detector  prevented  an 
experimental  validation  of  MCNP6.  The  semiconductor  properties,  and  not  the  material 
absorption  cross-sections,  drive  the  energy-dependence  of  the  germanium  detector. 
Therefore,  the  dose  deposition  provided  by  MCNP6  did  not  match  the  output  of  the 
germanium  detector.  The  experimentally-derived  GPF  was  1 .202,  and  the  output  spectra 
are  shown  in  Figure  23. 


Initial  Validation,  HPGe  initial  Validation,  HPGe 

Cs-1 37,  Free  Field  cs-1 37.  Shielded 


ZOO  KO  400  800  800  0  100  ?t>0  SW  400 

Energy  Deposited  (keV)  Energy  Deposited  (keV) 


Figure  23.  Experimental  spectra  from  a  HPGe  detector  using  an  identical  source  in  an  unshielded 
(left)  and  shielded  (right)  configuration. 


Initial  Validation,  HPGe 


Energy  Deposited  (keV) 


Initial  Validation,  HPGe 


Energy  Deposited  (keV) 


Figure  24.  MCNP6  of  the  HPGe  detector  in  the  free  field  (left)  and  shielded  (right)  spectra  using  the 
dose-deposition  tally.  Note  the  clear  Compton  edge  and  continuum. 
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A  multinuclide  source  allowed  for  an  energy-based  efficiency  calibration;  this 
calibration  compared  a  flux  of  known  photon  peak  activities  to  detected  peak  activities. 
This  process  inherently  included  the  influence  of  material  photon  absorption  cross 
sections.  Figure  24  shows  the  MCNP6  simulations’  dose-deposition  spectra.  Although  the 
662  keV  peak,  Compton  edge,  and  escape  peak  from  the  air  and  iron  were  at  the  correct 
energies,  the  shapes  of  the  models  differed.  The  proportions  of  the  modeled  features  in 
Figure  24  did  not  match  up  to  the  proportions  of  their  experimental  counterparts  in  Figure 
23.  This  resulted  in  a  GPF  of  1.042,  nearly  15%  less  than  the  experimental  results. 

A  logical  next  step  was  to  simply  tally  the  photon  flux  into  the  model  detector  and 
apply  an  experimentally-derived  efficiency  calibration.  This  solution  did  not  yield  usable 
results  either.  Although  the  photon  peak  was  adjusted  appropriately  to  account  for  the 
efficiency  calibration,  the  calibration  did  not  account  for  spectroscopic  features 
associated  with  Compton  scattering  inside  the  detector.  Without  a  Compton  edge  or 
continuum,  the  result  did  not  include  the  physics  necessary  to  compute  a  valid  dose  or 
GPF.  Figure  25  shows  the  modeled  flux  spectrum  with  the  energy  calibration  applied. 
Although  this  does  show  the  escape  peak  from  Compton  interactions  in  the  iron  and  air, 
the  lack  of  lower-energy  Compton  features  from  interactions  inside  the  detector  reduced 
the  resulting  GPF  to  0.977  for  this  approach,  23%  lower  than  the  experimental  results. 
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Simulated  Free  Field  Flux  into  HPGe  Simulated  Shielded  Flux  into  HPGe 

MCNP6  MCNP6 


Figure  25.  Simulated  photon  flux  spectra  from  the  free  field  (left)  and  shielded  (right)  configurations, 
adjusted  with  an  experimentally-derived  efficiency  calibration.  Note  the  lack  of  Compton  edges  or 
continua. 

The  fundamental  issue  with  validating  MCNP6  with  a  gennanium  detector  comes 
from  the  nature  of  the  semiconductor  detector  itself  Although  the  high  resolution  of  this 
detector  allows  for  precise  peak  identification,  its  efficiency  is  a  function  of  several 
factors,  including  the  photon  attenuation  properties  of  germanium,  the  electronic 
properties  of  the  germanium  for  moving  charges,  the  locations  and  sizes  of  the  p  and  n- 
type  regions  within  the  detector,  the  geometry  of  the  detector,  the  temperature  of  the 
Sterling  engine-cooled  detector,  and  the  stability  of  the  bias  power  supply. 

Figure  26  shows  an  experimentally-derived  energy-dependent  efficiency  of  the 
gennanium  detector.  This  efficiency  curve  includes  the  relevant  factors  and  can  be 
applied  to  a  known  source  or  flux  to  predict  the  output  of  the  detector.  However,  it 
yielded  invalid  results  when  applied  to  the  simulated  spectra. 
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HPGe  Detector  Efficiency 
Log-Log  Interpolation  of  Experimental  Data 


Figure  26.  Experimental  efficiency  calibration  for  the  germanium  detector. 


In  order  to  validate  MCNP6  for  dose  depositions  inside  a  germanium  detector,  the 
geometric  and  X-ray  attenuation  properties  would  have  to  be  segregated  into  a  discrete 
efficiency  function  from  the  electronic  efficiency  of  the  detector.  This  would  require  the 
development  of  a  signal  processing  algorithm,  and  that  algorithm  would  also  have  to  be 
verified  and  validated  experimentally.  While  such  a  study  could  be  a  useful  topic  for 
future  research,  it  was  beyond  the  scope  of  this  project.  Since  the  Nal  detector  offered  a 
more  straightforward  approach  for  validation,  it  was  the  only  detector  used  in  the 
remainder  of  the  observations  taken. 


5.5  Modeling  Resolution 

5.5.1  Geometry  Resolution 
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Since  the  initial  validation  was  successful,  a  series  of  alterations  reduced  the 
geometric  fidelity  of  the  model  to  the  experiment  to  find  a  “breaking  point.”  These 
alterations  were  grouped  together  into  two  categories.  Each  category  was  evaluated  at 
three  levels  of  fidelity. 

The  first  category  of  model  resolution  was  geometrical.  Since  the  initial  validation 
was  valid,  its  input  deck  was  adopted  as  the  “high  fidelity”  run.  The  “medium  fidelity” 
run  removed  the  far-side  plate  of  the  box  and  eliminated  all  air  outside  the  box.  The  “low 
fidelity”  eliminated  all  but  the  single  panel  of  the  box  between  the  source  and  the  Nal 
detector  and  replaced  all  air  with  vacuum. 

Table  3.  GPFs  produced  in  MCNP6  with  varying  degrees  of  geometric  fidelity  to  the  problem. 


Fidelity 

GPF 

±2  G 

High 

1.021 

0.022 

Medium 

1.031 

0.022 

Low 

1.098 

0.024 

The  GPFs  in  Table  1  produced  by  the  medium  and  low  resolution  models  differed 
from  that  produced  by  the  high  fidelity  model.  However,  when  compared  to  the 
experimental  GPF  of  1.0587  from  the  initial  validation  experiment,  all  three  models  met 
the  validation  criteria.  While  MCNP6  produced  virtually  identical  dose  deposition  spectra 
for  all  three  fidelity  levels  in  the  free  field  configurations,  the  shielded  spectra  of  Figure 
27  differ  significantly  in  the  low  fidelity  case.  This  suggests  that  a  reliable  GPF  may  be 
obtainable  from  geometric  approximations  of  an  actual  vehicle,  but  a  simple  barrier 
would  not  sufficiently  model  GPFs  of  a  three-dimensional  body,  such  as  an  annored 
vehicle. 
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MCNP6  High-Fidelity  Geometry 


Energy  Deposited  (keV)  per  Count 


MCNP6  Low-Fidelity  Geometry 
Shielded  Cs-1 37 
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Energy  Deposited  (keV)  per  Count 


Figure  27.  Comparison  of  MCNP6  spectra  with  geometries  in  the  high-fidelity  (left)  and  low  fidelity 
(right)  settings.  The  windows  of  both  are  cutting  off  the  full  height  of  the  full-energy  peaks  in  order 
to  highlight  the  differences  in  the  lower  energy  ranges. 


5. 5. 2  Material  Data  Resolution 

The  second  category  of  model  resolution  studied  focused  on  the  material  data 
cards.  The  “high  fidelity”  evaluation  included  each  element’s  natural  isotopes  in  their 
normal  relative  abundance.  The  “medium  fidelity”  evaluation  simplified  each  element  to 
its  most-common  isotope.  The  “low  fidelity”  evaluation  eliminated  all  but  the  most 
common  material  in  each  cell  with  the  exception  of  the  Nal,  which  contains  an  even  mix 
of  sodium  and  iodine  atoms. 
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Table  4.  GPFs  produced  in  MCNP6  with  varying  degrees  of  material  fidelity  to  the  problem. 


Fidelity 

GPF 

±2o 

High 

1.021 

0.016 

Medium 

1.021 

0.015 

Low 

1.023 

0.017 

Varying  the  material  data  fidelity  had  little  effect  on  the  GPF.  The  spectra 
produced  were  virtually  identical.  The  only  noticeable  spectral  difference  was  a 
sharpening  of  low-energy  escape  peaks  in  the  shielded  configuration  as  the  iron  was 
simplified  to  a  monoisotopic  approximation.  While  Figure  28  illustrates  this  difference,  it 
also  shows  little  change  to  the  overall  shape  of  the  spectrum.  It  follows  that  simplifying 
the  materials  has  little  effect  upon  the  overall  dose  or  GPF. 


MCNP6  High-Fidelity  Material  Data 


Energy  Deposited  (keV)  per  Count 


Figure  28.  MCNP  Shielded  Spectra  with  high  material  fidelity  (left)  and  low  material  fidelity  (right). 
The  shape  of  the  escape  peaks  around  200keV  is  the  only  significant  difference. 


5.5.3  Statistical  Resolution 

The  final  category  of  model  resolution  studied  the  influence  of  the  track  count  on 
the  model’s  GPF.  The  “high  fidelity”  evaluation  had  100  million  tracks  and  tallied  energy 
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depositions  per  track  in  1  keV  bins.  The  “medium  fidelity”  evaluation  had  10  million 
tracks  and  tallied  energy  depositions  per  track  in  10  keV  bins.  The  “low  fidelity” 
evaluation  had  100,000  tracks  and  tallied  energy  depositions  per  track  in  100  keV  bins. 

The  track  numbers  of  the  simulations,  the  run  times  of  the  experiments,  and  the 
bin  widths  of  both  drove  the  statistical  variance  of  the  MCNP6  models  and  experimental 
spectra.  MCNP6  output  files’  dose  deposition  tallies  provided  the  energy,  activity,  and 
relative  error  of  each  bin  in  its  output  file.  The  relative  error  is  defined  as  [22] 

cr- 

Re  lativeError  =  -=?- 
x 

where  a-  is  the  standard  deviation  for  that  bin’s  tally  and  xis  the  mean  activity  for  that 

bin.  The  relative  error  describes  the  precision  of  an  individual  energy  bin,  and  not  the 
overall  dose,  and  the  overall  dose  is  the  value  of  interest  to  this  work.  An  effective 
method  for  evaluating  the  overall  error  for  a  GPF  is  to  observe  the  highest  peak  bins  in 
the  free  field  and  shielded  simulations  of  a  GPF  calculation.  Multiplying  the  higher  of 
these  two  relative  errors  by  the  GPF  will  yield  a  worst-case  approximation  of  the  GPF’s 
standard  deviation. 

The  data  below  suggests  an  inverse  square  root  relationship  between  the  dose’s 
standard  deviation  and  the  track  count.  Since  95%  confidence  is  the  threshold  for 
validation  for  this  work,  the  modeled  GPF  ±  2  standard  deviations  should  lie  entirely 
within  the  measured  GPF  ±  5%  according  to  the  68-95-99.7  rule  [6].  The  initial 
validation’s  10  million  tracks  is  high  enough  to  meet  the  validation  standard,  while  a  run 
of  1  million  tracks  was  not,  as  its  confidence  interval’s  lower  bound  was  below  95%  of 
the  experimental  GPF. 
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Table  5.  Statistical  data  from  the  MCNP6  track  count  study.  Track  counts  meeting  the  validation 
criteria  must  have  their  confidence  intervals  bounded  between  the  experimental  GPF  ±  5%,  or 
between  1.0058  and  1.1116 


Number  of 
Tracks 

Free  Field  % 
Relative  Error 

Shielded  % 
Relative  Error 

GPF 

Lower 

Bound 

Confidence 

Interval 

Upper 

Bound 

Confidence 

Interval 

100000000 

0.29 

0.0034 

1.023 

1.016044 

1.029956 

10000000 

1.08 

0.009 

1.0203 

1.001935 

1.038665 

1000000 

0.0267 

0.034 

1.0053 

0.93694 

1.07366 

100000 

0.0924 

0.1066 

0.9413 

0.740615 

1.141985 

Although  the  initial  validation  experiment  was  successful,  and  closely  matched 
the  shape  of  the  experiment  with  an  appropriate  point  spread  function  applied,  the 
spectrum  produced  by  the  initial  MCNP6  simulation  had  bins  10  keV  wide,  which  was 
adequate  to  match  the  low-resolution  results  of  the  Nal  detector.  However,  MCNP6  can 
tally  in  bins  far  smaller  than  10  keV.  In  an  effort  to  achieve  a  closer  GPF  match  to 
experimental  results,  the  input  deck  was  adjusted  to  output  dose  depositions  in  1  keV 
increments,  while  the  number  of  tracks  was  increased  to  100  million.  With  these  settings, 
MCNP6  produced  a  GPF  of  1.022.  Since  this  result  matches  what  was  obtained  using  a 
lower-fidelity  approach,  this  higher-resolution  approach  was  of  little  use  for  calculating 
GPFs.  However,  this  high-resolution  technique  produced  spectral  features  that  confirmed 
MCNP6’s  ability  to  model  different  mechanisms  of  photon  transport,  which  are 
addressed  in  the  next  section. 
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5.6  Confirmation  of  Characteristic  X-rays  and  Pair  Production 


An  interesting  feature  appeared  when  evaluating  MCNP6  spectra  at  this  higher 
resolution  and  track  count.  When  using  a  dose  deposition  tally  of  1  keV  bin  increments 
for  100  million  tracks  or  more  in  a  cylinder  of  Nal  using  a  60Co  source,  four  anomalous, 
small  peaks  appear  in  the  resulting  spectrum  of  Figure  29.  Two  of  them  can  be  explained 
by  pair  production:  an  escape  peak  from  the  iron  box  at  511  keV,  and  an  escape  peak 
from  the  detector  at  823  keV,  which  is  511  keV  less  than  the  1332  peak.  However,  pair 
production  cannot  explain  why  a  small  but  distinct  a  “shadow  peak”  occurred  below  each 
known  peak.  In  every  case,  the  shadow  peak  was  over  an  order  of  magnitude  greater  in 
activity  than  the  bins  to  its  left  and  right,  and  it  was  precisely  28  keV  less  than  the  main 
peak. 


MCNP6  High-Resolution  Simulation 


Energy  Deposited  (keV)  per  Count 


Figure  29.  High-resolution  MCNP6  spectrum  for  6°Co  in  a  shielded  environment.  Note  the  “shadow 
peaks”  at  fixed  intervals  below  both  characteristic  peaks  for  the  6°Co. 
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The  “shadow  peak”  is  an  escape  peak  caused  by  a  characteristic  X-ray  escaping 
the  Nal.  When  a  peak  photon  deposits  all  its  energy  inside  the  Nal  detector,  it  can  excite 
an  L-shell  electron  of  the  iodine  into  the  K  shell.  This  process  emits  a  28  keV  X-ray  [5]. 
When  this  photon  escapes,  it  diminishes  the  resulting  peak  energy  deposited  by  28  keV. 
Although  the  escape  peak’s  low  activity  prevents  a  significant  effect  on  the  overall  GPF, 
its  presence  demonstrates  MCNP6’s  ability  to  model  material-dependent  characteristic  X- 
rays. 

This  higher-resolution  approach  also  allows  for  a  more  definitive  identification  of 
peaks  associated  with  pair  production  events.  Increasing  the  energy  of  the  photons  to 
4  MeV  drove  pair  production  to  become  a  dominant  effect  in  the  dose  deposition 
spectrum.  Figure  30  demonstrates  this  by  convoluting  the  MCNP6  spectrum  with  the  Nal 
detector  response.  Escape  peaks  at  51 1  and  1022  keV  less  than  the  primary  peak,  as  well 
as  the  511  keV  escape  peak,  demonstrate  pair  production’s  dominance  as  a  photon 
transport  mechanism  in  this  scenario. 
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MCNP6  High-Energy  Simulation 
4  MeV  Source,  Shielded 

xIO  1 

3.5- 


3 


Energy  Deposited  (keV)  per  Count  *  io‘ 

Figure  30.  MCNP6  Dose  Deposition  Spectrum  for  a  modified  version  of  the  Initial  Validation  model. 
The  primary  photon  energy  has  been  increased  to  4.0  MeV,  and  the  bins  are  1  keV  wide.  The 
MCNP6  output  has  been  convoluted  with  a  point  spread  function  to  simulate  an  expected  Nal  output. 

Practical  concerns  impeded  the  experimental  validation  of  these  transport 
mechanisms.  The  Nal  detector’s  large  size  prevented  most  characteristic  X-rays  from 
escaping,  and  its  resolution  was  too  low  to  differentiate  between  the  true  peak  and  the 
escape  peak.  Validating  pair  production  was  complicated,  as  the  only  4  MeV  source 
readily  available  at  AFIT  was  a  PuBe  source,  which  also  produces  significant  amounts  of 
neutrons  and  other  gamma  rays. 

Although  both  of  these  mechanisms  might  merit  further  investigation  for 
validating  photon  transport,  GPF  validation  required  validated  aggregate  results,  not 
specific  mechanisms.  Therefore,  neither  effect  merited  further  investigation  in  the  context 
of  this  work. 
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5. 7  Alternate  Source  Evaluations 


Following  the  initial  validation  experiment,  several  sources  were  substituted  for 
the  Cs  using  the  Nal  detector  and  the  iron  box  in  its  %  inch  wall  thickness 
configuration,  and  these  sources  were  simulated  in  MCNP6.  For  the  higher-energy 
sources,  the  MCNP6  GPFs  met  the  validation  criteria  when  compared  to  the  experimental 
GPFs. 

Table  6.  Alternate  source  GPF  Evaluations. 


Source 

High 

Peak 

Modeled 

GPF 

Model 
GPF 
±2  a 

Experimental 

GPF 

Experimental 
GPF  ±2  a 

60Co 

1332 

1.052 

0.0231 

1.054 

0.0116 

137Cs 

662 

keV 

1.022 

0.0225 

1.068 

0.0117 

(Jt 

O 

o 

123 

keV 

2.457 

0.0541 

4.065 

0.0447 

109Cd 

88  keV 

6.89 

0.152 

3.45 

0.0380 

The  57Co  and  109Cd  source  experiments  did  not  validate  their  respective  models. 
This  is  due  to  the  vast  majority  of  the  scattered  photons  produced  by  these  sources  falling 
under  80  keV,  which  was  too  low  for  the  Nal  detector’s  threshold. 

One  interesting  result  from  this  study  was  that  60Co  GPF  approximately  matched 
(and  slightly  exceeded  in  the  model)  the  Cs  GPF,  despite  its  much-higher  photon 
energy  and  corresponding  lower  material  cross  sections.  The  energy  dependence  of  the 
Klein-Nishina  distribution  explains  this  behavior.  At  the  1173  and  1332  keV  energies 
emitted  by  60Co,  Compton  scattering  will  tend  to  deflect  photons  at  small  angles,  while 
the  lower  energies  of  Cs  tends  to  deflect  at  wider  angles.  The  lower-energy  photons  of 
the  ~  Cs  will  tend  to  scatter  off  the  walls  of  the  walls  of  the  box  toward  the  detector  more 
often  than  the  primary  photons  coming  from  the  60Co. 
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5.8  Thick-Wall  Assembly  Experiments  and  Models 


The  iron  box  was  then  reassembled  to  its  maximum  wall  thickness  of  1%  inches 
and  tested  in  a  variety  of  configurations  to  demonstrate  MCNP6’s  ability  to  model 
scenarios  with  complicated  geometries  and  sources.  Due  to  the  need  for  a  stronger  table 
for  the  box  and  better  shielding  for  the  longer  observation  times  required,  the  box  was 
moved  into  a  smaller,  concrete-enclosed  room.  This  caused  significantly  more  lower- 
energy  activity  in  the  “free  field”  measurements  due  to  scattering  off  the  walls,  so  the 
room  was  modeled  as  a  concrete  box  around  a  box  of  air  containing  the  steel  enclosure. 
Since  little  of  the  lower-energy  sources’  photons  penetrated  the  thin- wall  configuration, 
only  137Cs  and  60Co  were  used  for  evaluating  the  thick-wall  assembly.  Sources  were 
placed  either  outside  the  center  of  a  face,  outside  the  center  of  the  edge,  or  outside  a 
comer,  as  shown  in  Figure  3 1 . 


Figure  31.  Placement  of  sources  for  the  thick-assembly  observations  and  models.  FC  is  face-centered, 
EC  is  edge-centered,  and  CN  is  corner.  The  source  is  placed  14  inches  from  the  center  of  the  detector 
in  along  one,  two,  and  all  three  axes  of  the  box,  respectively. 

Figure  32  and  Figure  33  demonstrate  how  the  experimental  spectra  for  the  free 
field  and  shielded  configurations  of  the  thick  box  experiments  closely  matched  their 
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MCNP6  simulations.  See  Appendix  A  for  the  other  thick-wall  assembly  experimental  and 
MCNP6  spectra.  Similar  to  the  initial  validation  experiment,  the  discrepancies  were 
most-noticeable  in  the  lower-energy  regions  of  the  free  field  observations,  where  the 
experimental  spectra  showed  more  activity  than  the  models.  As  noted  in  Section  5.3,  this 
is  likely  driven  by  characteristics  of  the  detector  system. 

The  first  evaluation  was  a  series  of  observations  similar  to  the  initial  validation 
experiment.  This  evaluation  produced  GPFs  for  a  1  l7Cs  and  60Co  source  just  outside  the 
center  of  the  box’s  face.  Table  7  contains  the  experimental  GPFs  and  their  corresponding 
GPFs  produced  in  MCNP6. 


Cs-137,  Face-Centered  Placement 
Free  Field  MCNP6  Spectrum  (Thick  Box  Series) 


Energy  Deposited  (keV) 


Cs-1 37,  Face-Centered  Placement 
Free  Field  Experimental  Spectrum  (Big  Box  Series) 


Energy  Deposited  (keV) 


Figure  32.  Comparison  of  free  field  MCNP6  simulation  (left)  and  the  corresponding  experimental 

137 

observation  (right),  for  the  thick-box,  face-centered  Cs  evaluation. 


Cs-137,  Face-Centered  Placement 
Shielded  MCNP6  Spectrum  (Thick  Box  Series) 


Energy  Deposited  (keV) 


Cs-137,  Face-Centered  Placement 
Shielded  Experimental  Spectrum  (Big  Box  Series) 


Energy  Deposited  (keV) 


Figure  33.  Comparison  of  shielded  MCNP6  simulation  (left)  and  the  corresponding  experimental 

137 

observation  (right),  for  the  thick-box,  face-centered  Cs  evaluation. 
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Table  7.  Experimental  and  MCNP  Gamma  Protection  Factors  for  face-centered  sources  outside  the 
thick-assembly  box. 


Source 

Os 

O 

O 

o 

137Cs 

Model  GPF 

1.76 

2.16 

Model  GPF  ±2a 

0.034 

0.042 

Experimental  GPF 

1.83 

2.12 

Experimental  ±2o 

0.013 

0.018 

The  thicker  walls  attenuate  the  peak-energy  photons  more  and  block  much  of  the 
down-scattered  low-energy  build-up  region.  While  this  lower-energy  region  remains 
active  in  the  shielded  configurations,  the  thicker  iron  results  in  fewer  photons  escaping 
from  the  iron  after  undergoing  Compton  interactions.  Thus,  the  thicker  walls  result  in 
significantly  higher  GPFs  than  in  the  thin-wall  assembly  for  both  the  60Co  and  the  137Cs 
as  shown  in  Table  7.  By  making  the  walls  five  times  thicker,  the  GPF  for  these  higher- 
energy  photons  has  gone  from  a  nearly  worthless  protection  factor  very  close  to  1  to  a 
better,  if  still  modest  value. 

The  second  evaluation  challenged  MCNP6  to  model  a  single  source  placed 
outside  an  edge  of  the  box.  This  increased  the  straight-line  distance  required  for  an 
attenuated  primary  photon  to  penetrate  the  steel  box.  It  also  increased  the  factors 
associated  with  photons  striking  a  face  of  the  iron  and  scattering  through  it  towards  the 
detector.  Observations  were  taken  in  a  free  field  and  shielded  configuration,  and  the 
resulting  GPFs  are  shown  in  Table  8. 
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Table  8.  Experimental  and  MCNP  Gamma  Protection  Factors  for  edge-centered  sources  outside  the 
thick-assembly  box. 


Source 

ON 

O 

n 

o 

137Cs 

Model  GPF 

2.89 

4.47 

Model  GPF  ±2cr 

0.055 

0.086 

Experimental  GPF 

3.29 

5.11 

Experimental  ±2a 

0.018 

0.028 

1  T7 

This  “edge  in”  configuration  significantly  increased  the  GPFs  for  both  Cs  and 
60Co.  Although  the  straight-line  path  through  the  iron  was  only  40%  greater  than  the  face- 
centered  configuration,  the  GPFs  were  approximately  twice  as  high.  This  is  likely  due  to 
the  geometric  factors  associated  with  the  setup.  While  nearly  50%  of  the  photons  from  a 
face-centered  source  will  strike  the  box,  only  25%  of  an  edge-centered  source  will  strike 
the  box.  Additionally,  all  of  the  photons  that  do  strike  the  panels  will  do  so  at  or  below  an 
incident  angle  of  45  degrees,  which  increases  their  path  length  in  the  iron.  The 
experimental  GPFs  were  both  higher  than  the  MCNP6  models  by  12%.  It  is  likely  that  the 
aluminum  frame,  which  lay  on  the  straight-line  path  between  the  source  and  the  detector, 
is  a  strong  contributor  to  this  discrepancy.  The  Mathematica  peak  attenuation  model 
demonstrates  that  the  aluminum  frame  should  contribute  a  3.01%  GPF  gain  for  60Co  and 
a  4.27%  GPF  gain  for  ~  Cs  peak  photons.  The  shorter  mean  free  paths  of  downscattered 
photons  from  the  steel  would  likely  drive  this  GPF  up  further.  However,  due  to 
complications  with  the  model  geometry  inputs,  the  frame  was  not  simulated;  this  would 
be  a  good  area  for  continued  research. 

The  third  evaluation  placed  the  source  outside  a  “point  comer”  of  the  box  in  an 
effort  to  further  challenge  its  ability  to  handle  three-dimensional  photon  transport.  By 
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placing  the  source  along  a  straight  line  from  the  center  of  the  detector  in  the  center  of  the 
box  through  a  corner  of  the  box,  this  evaluation  maximized  the  straight-line  distance 
required  for  an  attenuated  primary  photon  to  penetrate  the  steel  box  and  interact  inside 
the  detector.  It  also  maximized  the  factors  associated  with  photons  striking  a  face  of  the 
iron  and  scattering  through  it  towards  the  detector. 


Table  9.  Experimental  and  MCNP  Gamma  Protection  Factors  for  sources  outside  the  corner  of  the 
thick-assembly  box. 


Source 

On 

O 

O 

o 

137Cs 

Model  GPF 

4.28 

8.96 

Model  GPF  ±2c> 

0.085 

0.177 

Experimental  GPF 

4.88 

9.60 

Experimental  GPF 
±2o 

0.043 

0.090 

The  corner  approach  offered  the  best  GPFs  of  the  three  configurations,  shown  in 
Table  9  to  be  approximately  three  times  more  than  the  face-centered  configuration.  This 
is  due  to  geometric  factors  which  reduce  the  number  of  photons  incident  upon  the  box. 
While  almost  half  of  the  emitted  photons  in  the  face-centered  configuration  strike  the 
box,  slightly  less  than  an  eighth  strike  it  in  the  comer-in  configuration.  Thus,  fewer 
primary  photons  have  a  chance  to  scatter  towards  the  detector,  and  most  of  those  that  do 
must  pass  through  an  increased  thickness  of  iron  due  to  their  oblique  incident  angle  to  the 
panels  of  the  iron  box. 

Again,  the  aluminum  frame  lay  in  the  direct  path  between  the  source  and  the 
detector,  so  the  experimental  GPF  was  notably  higher  than  the  MCNP6  model. 
Interestingly,  the  degree  to  which  the  experimental  GPF  exceeded  the  model  was  by  a 
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consistent  0.7.  However,  the  proportional  difference  was  12%  for  the  60Co  and  6%  for  the 
Cs,  since  the  Cs  GPF  was  so  significantly  higher.  Part  of  these  discrepancies  could 
be  accounted  for  by  the  extra  peak  attenuation  of  the  aluminum,  which  increased  the 
attenuation  by  5.27%  for  137Cs  and  3.70%  for  60Co.  However,  the  downscattered  photons 
required  the  aluminum  frame  to  be  modeled  in  MCNP6,  and  the  input  deck  was  difficult 
to  visualize  with  available  software. 

The  final  evaluation  was  designed  to  test  MCNP6’s  ability  to  model  multiple 
sources  from  multiple  directions  simultaneously.  A  60Co  source  was  placed  outside  the 

i  'in 

comer  of  the  box,  while  a  Cs  source  was  placed  at  the  center  of  an  edge.  Observations 
were  taken  in  a  free  field  and  shielded  configuration. 


Multi-Source,  Cs-137  at  Center-Edge,  Co-60  at  Corner  Multi-Source,  Cs-137  at  Center-Edge,  Co-60  at  Comer 

.  Free  Field  MCNP6  Spectrum  (Thick  Box  Series)  Free  Field  Experimental  Spectrum  (Thick  Box  Series) 


Energy  Deposited  (keV)  Energy  Deposited  (keV) 


Figure  34.  Comparison  of  free  field  MCNP6  simulation  (left)  and  the  corresponding  experimental 
observation  (right),  for  the  thick-box,  multi-source  evaluation. 


Multi-Source,  Cs-1 37  at  Center-Edge,  Co-60  at  Corner  Multi-Source,  Cs-1 37  at  Center-Edge,  Co-60  at  Corner 

Shielded  MCNP6  SDectrum  (Thick  Box  Series!  Shielded  Experimental  Spectrum  (Thick  Box  Series) 


5.M  tCOO  1500  0  KB  1CCO 

Energy  Deposited  (keV)  Energy  Deposited  (keV) 


Figure  35.  Comparison  of  shielded  MCNP6  simulation  (left)  and  the  corresponding  experimental 
observation  (right),  for  the  thick-box,  multi-source  evaluation. 
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Table  10.  Gamma  Protection  Factors  for  Multi-Source  Experiment  and  Model. 


GPF 

2sigma 

Experiment 

5.03 

0.031 

Model 

4.6 

0.101 

Using  available  reference  data  for  the  radioisotope  sources,  the  current  activities 
of  the  137Cs  and  60Co  sources  were  calculated;  these  activities  were  then  used  to  provide 
the  appropriate  probabilities  of  each  peak  in  MCNP6.  The  resulting  spectra  were  a  very 
close  match  to  MCNP6,  as  demonstrated  in  Figure  34  and  Figure  35.  Table  10  shows  that 
the  experimental  GPF  slightly  exceeds  the  value  of  the  model,  most  likely  due  to 
additional  protection  provided  by  the  un-modeled  aluminum  frame. 

5. 9  Sources  of  Modeling  Error 

Although  Section  5.9  addresses  the  statistical  error  for  the  work,  there  were  many 
error  sources  that  were  difficult  to  quantify  for  both  the  model  and  the  experiment.  Most 
of  these  sources  were  assumed  to  have  little  impact. 

While  Section  5.5  explores  the  role  of  model  fidelity,  it  does  not  evaluate  a 
number  of  modeling  assumptions  and  approximations.  All  MCNP6  models  treated  the 
sources  as  point  sources,  rather  than  as  the  round  planchettes  they  are.  The  geometry  of 
the  room  was  either  not  included  for  thin  assembly  simulations  or  simplified  to  be  a 
simple  cube  for  the  thick  assembly  simulations.  The  model  included  none  of  the  other 
equipment  in  the  room,  including  the  table  the  box  was  on. 
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The  steel  box  also  had  many  simplifications.  The  box’s  aluminum  frame  and  bolts 
were  not  included  in  the  model,  nor  were  the  holes  in  the  plates  for  the  bolts  and  the 
detector  cables.  This  likely  drove  the  experimental  GPF  slightly  higher,  especially  for 
experiments  with  sources  outside  the  edges  and  comers.  The  thin  assembly’s  plates 
overlapped  at  the  corners,  but  this  overlap  was  not  modeled,  since  sources  for  this 
assembly  were  always  placed  at  the  center  of  a  plate’s  face.  The  thick  assembly’s  plates 
consisted  of  multiple  powder-coated  plates,  and  these  coatings  were  not  modeled,  as  they 
were  assumed  to  be  of  insufficient  mass  to  be  significant  contributors  to  gamma 
attenuation. 

The  detector  was  modeled  as  a  single  cylinder  of  Nal.  The  model  did  not  include 
its  aluminum  sheath,  photo-multiplier  tube,  or  associated  cables.  The  model  also  left  out 
the  aluminum  stand  which  held  the  detector  in  the  middle  of  the  iron  box. 

The  fidelity  investigation  of  Section  5.5  suggests  that  these  approximations  likely 
had  little  effect  on  the  overall  work.  If  a  back  wall  of  the  box  can  be  removed  with 
minimal  effect,  it  follows  that  small  holes  or  extra  volumes  made  of  low-mass  gasses 
should  also  have  a  minimal  effect. 

5.10  Sources  of  Experimental  Error 

Low-level  noise  in  the  detector  caused  the  greatest  amount  of  experimental  error. 
The  upper  limit  for  this  noise  varied  depending  on  the  source  and  shielding  between  55  to 
85  keV,  and  the  low-level  discriminator  had  to  be  adjusted  to  avoid  having  this  noise 
flood  the  analog-digital  converter  with  meaningless  data.  The  likely  cause  of  this  low- 
level  noise  is  fluorescence  in  the  Nal  which  tends  to  build  up  and  release  over  time  from 
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the  scintillator.  As  a  result,  none  of  the  photons  from  0  to  approximately  80  keV  were 
counted  in  the  experiments.  Since  these  photons  tended  to  be  shielded  more  effectively 
than  the  higher  energy  photons,  this  tended  to  drive  the  experimental  GPFs  down. 
However,  since  most  of  the  actual  dose  came  from  higher-energy  photons,  this  had  a 
relatively  small  effect  on  the  overall  dose  calculations  and  GPFs. 

The  high-voltage  bias  tended  to  fluctuate  by  approximately  20  V.  Although  the 
detector  response  function  spread  out  peaks  far  more  than  these  voltage  fluctuations 
would  do  on  their  own,  the  fact  remains  that  this  likely  contributed  to  the  spreading. 
Future  work  could  be  improved  by  using  a  power  regulation  device. 

When  fully  assembled  with  1  %  inch  walls,  the  bolts  holding  the  box  together  were 
not  a  tight  fit,  nor  were  there  enough  bolts  to  fill  every  hole  in  every  plate.  The  resulting 
mechanical  play  caused  offsets  of  a  few  millimeters  and  gaps  of  up  to  a  millimeter 
between  adjacent  plates.  These  imperfections  had  little  influence  on  the  overall 
experiment,  but  they  could  be  a  factor  for  future  experiments  with  this  box,  particularly  if 
a  source  is  positioned  so  that  it  partly  “shines  through”  one  of  these  holes. 

Finally,  the  dose  calculations  for  the  experimental  measurements  assumed  a  linear 
relationship  between  energy  and  bins  for  the  Nal  detector.  Since  the  dose  preferentially 
weights  high-energy  peaks,  the  highest-energy  peak  in  each  spectra  detennined  the 
calibration  of  the  entire  spectrum.  Although  the  detector’s  response  is  close  to  linear, 
these  calculations  could  be  more  precise  if  they  applied  an  interpolation  between  known 
calibration  peaks. 
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6.  Conclusions 


6.1  Verifying  MCNP6  Against  Other  Options 

MCNP6  compares  favorably  to  simplified-geometry  simulations  and  GEANT4. 
The  Mathematica  model,  the  MATLAB  Monte  Carlo  script,  and  GEANT4  were 
markedly  slower  and  less  portable.  The  custom  Mathematica  and  MATLAB  models 
could  not  approach  the  comprehensive  particle  simulation  capabilities  of  MCNP6.  While 
GEANT4  offered  comparable  capabilities,  its  higher  input  file  complexity  and  slower 
computational  time  made  it  less  suitable  than  MCNP6,  and  the  spectrum  it  produced  did 
not  match  the  experimental  results  as  closely  as  the  spectrum  produced  by  MCNP6. 
While  other  particle  simulators  might  also  suffice  for  this  work,  the  comparisons 
perfonned  in  this  research  verify  MCNP6  as  an  excellent  platform  for  evaluating  GPFs. 

6.2  Validating  MCNP6  to  Produce  Gamma  Protection  Factors 

The  models  and  experiments  support  MCNP6’s  validation  for  producing  reliable 
GPFs.  The  program  was  able  to  produce  both  GPFs  and  spectra  that  closely  matched 
experimental  results  despite  significant  simplifications  and  approximations  in  the  model. 
The  limits  of  this  validation  were  not  driven  by  MCNP6  or  the  available  computing 
power,  but  rather  the  limitations  of  the  detector  systems  used.  The  Germanium  detector, 
while  offering  excellent  spectral  resolution,  was  inadequate  for  evaluating  GPFs  due 
energy  efficiency  factors  unrelated  to  its  photon-based  material  cross  sections.  The  Nal 
detector  was  a  far  superior  tool  for  GPF  validation,  though  its  performance  below  75  keV 
reduced  its  utility  for  evaluating  low-energy  sources. 
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6.3  Examining  Different  Influences  on  Gamma  Protection  Factors 


Varying  the  photon  energies  by  using  different  radioisotopes  showed  a  tendency 
towards  high  GPFs  below  150  keV,  even  with  relatively  thin  shielding.  GPFs  for  photon 
sources  at  662,  1173,  and  1332  suggest  a  relatively  flat  region  which  corresponds  to  the 
dominance  of  the  Compton  scatter  interaction.  In  this  region,  there  is  a  GPF-neutral 
tradeoff  between  the  lower  attenuation  of  higher-energy  photons  and  the  tendency  to 
scatter  at  wider  angles  for  lower-energy  photons. 

The  edge  and  comer  studies  in  the  thick-wall  experiments  suggest  that  shielding 
corners  provide  nearly  as  much  protection  as  a  flat  plate  with  the  same  thickness  as  the 
path  length  through  the  corner.  This  is  of  particular  interest  for  vehicles  fielded  during 
Operation  Iraqi  Freedom.  The  “V-shape”  bottom  of  their  hulls,  created  to  deflect  blast, 
should  also  offer  improved  protection  from  surface-deposited  radioactive  contamination 
directly  underneath  the  vehicle.  Similarly,  this  suggests  that  the  sloped  frontal  annor  of 
the  Abrams  tank  should  provide  far  better  protection  than  its  sides  or  back  against  the 
photon  radiation  burst  of  a  nuclear  weapon.  In  a  civilian  setting,  a  medical  photon-based 
radiation  treatment  center  might  benefit  from  having  its  operator  work  from  inside  a 
shielded  enclosure  with  the  edge  of  the  enclosure  pointing  toward  the  radiation  source. 

6.4  Future  Work 

The  series  of  experiments  on  the  thick-wall  assembly  suggest  a  relationship 
between  wall  thickness,  peak  energy,  geometric  configuration,  and  gamma  protection 
factors  that  warrants  further  investigation.  While  the  73%  increase  in  straight-line 
shielding  associated  with  moving  a  source  from  the  face  of  the  box  to  the  comer  roughly 
corresponds  to  a  proportional  increase  in  its  GPF,  the  proportionality  varies  with  energy. 
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Characterizing  this  relationship  would  require  a  more  extensive  range  of  photon  energies, 
and  could  suggest  better  shielding  geometries  against  known  sources  and  locations. 

Another  aspect  of  this  problem  that  could  produce  some  interesting  results  is 
evaluating  the  role  of  the  detector’s  position  within  the  annored  enclosure.  This  would 
likely  play  a  more  significant  role  for  penetrating  photons  that  are  prone  to  scattering  at 
wider  angles. 

Although  some  of  the  models  presented  confirmed  that  MCNP6  uses 
characteristic  X-rays  and  pair  production  in  its  photon  transport,  these  mechanisms  were 
not  experimentally  validated.  Although  the  experimental  spectra  can  validate  the 
program’s  ability  to  model  the  Compton  effect,  these  other  mechanisms  would  require 
different  sources  and  detectors. 

The  most  important  next  step  is  to  combine  the  work  of  this  research  with 
previous  work  in  validating  MCNP6  for  neutron  protection  factors  into  a  single  radiation 
protection  factor  (RPF)  validation  in  a  single  experiment  and  a  single  model.  This  could 
be  evaluated  using  the  fast  burst  reactor  at  White  Sands  Missile  Range  or  a  PuBe  source 
at  the  Air  Force  Institute  of  Technology. 
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Appendix  A:  Experimental  and  MCNP6  Spectra 

This  appendix  provides  the  experimental  and  modeled  spectra  taken.  In  most  cases,  the 
results  are  presented  as  side-by-side  comparisons  of  experimental  spectra  and  modeled 
spectra.  The  left  side  is  the  MCNP6  model  spectrum,  and  the  right  side  is  the 
experimental  spectrum. 

Each  figure  is  scaled  to  allow  for  a  shape  comparison  between  the  model  and 
experimental  spectra.  Unless  otherwise  noted,  the  modeled  spectra  have  been  convoluted 
with  a  point  spread  function  to  approximate  the  output  of  the  detector  system.  All  spectra 
were  taken  with  the  detector  centered  inside  the  steel  enclosure. 

Some  results  intended  to  study  the  characteristics  of  the  MCNP6  model  do  not  have 
corresponding  experimental  results. 


Initial  Validation 

1  T7 

Initial  validation  spectra  were  taken  with  a  Cs  source  14  inches  away  from  the  center 
of  the  detector.  In  the  shielded  configuration,  the  source  was  placed  outside  the  center  of 
one  of  the  enclosure’s  steel  plates. 


86 


Initial  Validation,  Free  Field  MCNP  Simulation,  10  Million  Tracks 
...  Cs-137  Source  15  Inches  from  Nal  Detector 


Initial  Validation,  Free  Field  Experiment 
Cs-1 37  Source  1 5  Inches  from  Nal  Detector 


Energy  of  Interactions  (keV) 


Figure  36.  Initial  validation,  comparison  of  free  field  MCNP6  simulation  (left)  and  the 
corresponding  experimental  observation  (right). 


Initial  Validation,  Shielded  MCNP6  Simulation,  10  Million  Tracks  Initial  Validation,  Shielded  Experiment 

Cs-137  Source  15  Inches  from  Nal  Detector  Cs-137  Source  15  Inches  from  Nal  Detector 


Energy  of  Interaction  (keV)  Energy  of  Interactions  (keV) 


Figure  37.  Initial  validation,  Initial  validation  with  1/4”  steel  walls,  comparison  of  MCNP6  simulation 
(left)  and  the  corresponding  experimental  observation  (right). 


Thin-Box  Observations  with  Different  Sources 

The  remainder  of  thin-box  evaluations  were  set  up  identically  to  the  initial 
validation  experiment  except  for  the  sources.  Different  sources  were  used,  though  their 
placement  was  the  same  as  the  initial  validation  experiment,  14  inches  from  the  detector, 
centered  on  the  face  of  one  of  the  enclosure’s  walls. 
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Figure  38. 
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Co  source,  comparison  of  free  field  MCNP6  simulation  (left)  and  the  corresponding 


experimental  observation  (right). 


MCNP6  Dose  Deposition  Spectrum 


Enerav  Deoosited  (keVl 


Dose  Deposition  Spectrum 
Shielded  Nal  Detector.  Cs-137 


Energy  Deposited  (keV) 


Figure  39. 
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Co  source,  comparison 


of  shielded  MCNP6 


simulation  (left)  and  the  corresponding 


experimental  observation  (right). 


MCNP6  Dose  Deposition  Spectrum 


Dose  Deposition  Spectrum 
Free  Field  Nal  Detector.  Co-57 


Energy  Deposited  (keV)  per  Count 
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Figure  40. '  Co  source,  comparison  of  free  field  MCNP6  simulation  (left)  and  the  corresponding 
experimental  observation  (right). 


MCNP6  Dose  Deposition  Spectrum 


Dose  Deposition  Spectrum 
Shielded  Nal  Detector,  Co-57 
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Energy  Deposited  (keV)  per  Count 
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Figure  41. ~  Co  source,  comparison  of  shielded  MCNP6  simulation  (left)  and  the  corresponding 
experimental  observation  (right). 


MCNP6  Dose  Deposition  Spectrum  Dose  Deposition  Spectrum 

Free  Field  Nal  Detector,  Cd-1 09  Free  Field  Nal  Detector,  Cd-1 09 


Figure  42.  109Cd  source,  comparison  of  free  field  MCNP6  simulation  (left)  and  the  corresponding 
experimental  observation  (right). 


MCNP6  Dose  Deposition  Spectrum  Dose  Deposition  Spectrum 

Free  Field  Nal  Detector,  Cd-1 09  Free  Field  Nal  Detector,  Cd-1 09 


Figure  43.  109Cd  source,  comparison  of  shielded  MCNP6  simulation  (left)  and  the  corresponding 
experimental  observation  (right). 


89 


Fidelity  Studies 


Fidelity  Studies  were  conducted  in  MCNP6  only  and  were  not  convoluted  with 
any  point  spread  functions. 


MCNP6  Hiah-Fidelitv  Geometry  MCNP6  Low-Fidelity  Geometry 


Energy  Deposited  (keV)  per  Count  Energy  Deposited  (keV)  per  Count 


Figure  44.  Comparison  of  high  and  low  geometric  fidelity  spectra. 


MCNP6  High-Fidelity  Material  Data 
Shielded  Cs-137 


Energy  Deposited  (keV)  per  Count 


MCNP6  Low-Fidelity  Material  Data 


Energy  Deposited  (keV)  per  Count 


Figure  45.  Comparison  of  high  and  low  material  data  fidelity  spectra. 
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Thick  Box  Observations 


These  observations  were  taken  with  the  box  in  its  full-thickness  1.25” 
configuration  using  137Cs  and  60Co  sources.  The  first  series  of  observations  and  models 


placed  the  source  outside  the  center  of  one  of  the  box’s  faces. 


Cs-137.  Face-Centered  Placement 


Free  Field  MCNP6  Spectrum  (Thick  Box  Series) 


Energy  Deposited  (keV) 


Cs-137,  Face-Centered  Placement 
Free  Field  Experimental  Spectrum  (Big  Box  Series) 


Energy  Deposited  (keV) 


Figure  46.  Comparison  of  free  field  MCNP6  simulation  (left)  and  the  corresponding  experimental 

137 

observation  (right),  for  the  thick-box,  face-centered  Cs  evaluation. 


Cs-137,  Face-Centered  Placement 
Shielded  MCNP6  Spectrum  (Thick  Box  Series) 


Energy  Deposited  (keV) 


Cs-137,  Face-Centered  Placement 
Shielded  Experimental  Spectrum  (Big  Box  Series) 

1Af - - - - — ^ - - 


Energy  Deposited  (keV) 


Figure  47.  Comparison  of  shielded  MCNP6  simulation  (left)  and  the  corresponding  experimental 

137 

observation  (right),  for  the  thick-box,  face-centered  Cs  evaluation. 
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Co-60,  Face-Centered  Placement  Co-60,  Face-Centered  Placement 

.  Free  Field  MCNP6  Spectrum  (Thick  Box  Series)  Free  Field  Experimental  Spectrum  (Thick  Box  Series) 


Energy  Deposited  (keV)  Energy  Deposited  (keV) 


Figure  48.  Comparison  of  free  field  MCNP6  simulation  (left)  and  the  corresponding  experimental 
observation  (right),  for  the  thick-box,  face-centered  6°Co  evaluation. 


Co-60,  Face-Centered  Placement  Co-60,  Face-Centered  Placement 

Shielded  MCNP6  SDectrum  (Thick  Box  Series'!  Shielded  MCNP6  Spectrum  (Thick  Box  Series) 


Energy  Deposited  (keV)  Energy  Deposited  (keV) 


Figure  49.  Comparison  of  shielded  MCNP6  simulation  (left)  and  the  corresponding  experimental 
observation  (right),  for  the  thick-box,  face-centered  6°Co  evaluation. 


Edge-Centered  Thick  Box  Observations 

The  second  series  of  observations  placed  the  source  outside  the  edge.  The  sources 
were  14  inches  away  along  the  axis  of  the  box’s  horizontal  edges. 


Cs-1 37,  Edge-Centered  Placement  Cs-1 37,  Edge-Centered  Placement 


Energy  Deposited  (keV)  Energy  Deposited  (keV) 


Figure  50.  Comparison  of  free  field  MCNP6  simulation  (left)  and  the  corresponding  experimental 
observation  (right),  for  the  thick-box,  edge-centered  6°Co  evaluation. 
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Cs-137,  Edge-Centered  Placement  Cs-137,  Edge-Centered  Placement 

Shielded  MCNP6  Spectrum  (Thick  Box  Series)  Free  Field  Experimental  Spectrum  (Thick  Box  Series) 


200  a»  400  S00  002  700  BOO  o  ICO  200  300  -CD  000 

Energy  Deposited  (keV)  Energy  Deposited  (keV) 


Figure  51.  Comparison  of  shielded  MCNP6  simulation  (left)  and  the  corresponding  experimental 

137 

observation  (right),  for  the  thick-box,  edge-centered  Cs  evaluation. 


Co-60,  Edge-Centered  Placement  Co-60,  Edge-Centered  Placement 


SCO  1000  15CO  0  500  ICC© 

Energy  Deposition  (keV)  Energy  Deposited  (keV) 


Figure  52.  Comparison  of  free  field  MCNP6  simulation  (left)  and  the  corresponding  experimental 
observation  (right),  for  the  thick-box,  edge-centered  6°Co  evaluation. 


Co-60,  Edge-Centered  Placement  Co-60,  Edge-Centered  Placement 


Energy  Deposited  (keV)  Energy  Deposited  (keV) 


Figure  53.  Comparison  of  shielded  MCNP6  simulation  (left)  and  the  corresponding  experimental 
observation  (right),  for  the  thick-box,  edge-centered  6°Co  evaluation. 
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Corner  Thick  Box  Observations 


The  third  series  of  observations  placed  the  source  outside  a  comer  of  the  box.  The 
sources  were  placed  14  inches  from  the  box’s  center  along  all  three  axes  of  the  box’s 
edges. 


Cs-1 37,  Corner  Placement  Cs-1 37.  Comer  Placement 


Energy  Deposited  (keV)  Energy  Deposited  (keV) 

Figure  54.  Comparison  of  free  field  MCNP6  simulation  (left)  and  the  corresponding  experimental 

137 

observation  (right),  for  the  thick-box,  corner  Cs  evaluation. 


Cs-1 37,  Comer  Placement 
,  Shielded  MCNP6  Spectrum  (Thick  Box  Series) 


Energy  Deposited  (keV) 


Cs-1 37,  Comer  Placement 
Shielded  Experimental  Spectrum  (Thick  Box  Series) 


Energy  Deposited  (keV) 


Figure  55.  Comparison  of  shielded  MCNP6  simulation  (left)  and  the  corresponding  experimental 


observation  (right),  for  the  thick-box,  corner 


137Cs 


evaluation. 
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Co-60,  Comer  Placement  Co-60,  Comer  Placement 


Energy  Deposited  (keV)  Energy  Deposited  (keV) 


Figure  56.  Comparison  of  free  field  MCNP6  simulation  (left)  and  the  corresponding  experimental 

137 

observation  (right),  for  the  thick-box,  corner  Cs  evaluation. 


Co-60,  Comer  Placement 


Energy  Deposited  (keV) 


Co-60,  Corner  Placement 

Shielded  Experimental  Spectrum  (Thick  Box  Series) 


Energy  Deposited  (keV) 


Figure  57.  Comparison  of  shielded  MCNP6  simulation  (left)  and  the  corresponding  experimental 


observation  (right),  for  the  thick-box,  corner 


60Co 


evaluation. 


Multi-Source  Observations 

The  final  series  of  observations  and  models  placed  two  sources  in  different  locations 
outside  the  enclosure.  The  enclosure’s  walls  for  this  observations  were  1.25  inches  thick. 
The  1 ,7Cs  was  placed  along  the  center  of  one  edge,  while  the  60Co  outside  a  comer.  These 
locations  were  identical  to  the  earlier  observations. 
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Multi-Source,  Cs-137  at  Center-Edge,  Co-60  at  Corner  Multi-Source,  Cs-137  at  Center-Edge,  Co-60  at  Comer 

.  Free  Field  MCNP6  Spectrum  (Thick  Box  Series)  Free  Field  Experimental  Spectrum  (Thick  Box  Series) 


Energy  Deposited  (keV)  Energy  Deposited  (keV) 


Figure  58.  Comparison  of  free  field  MCNP6  simulation  (left)  and  the  corresponding  experimental 
observation  (right),  for  the  thick-box,  multi-source  evaluation. 


Multi-Source,  Cs-137  at  Center-Edge.  Co-60  at  Corner 
Shielded  Experimental  Spectrum  (Thick  Box  Series) 


Figure  59.  Comparison  of  shielded  MCNP6  simulation  (left)  and  the  corresponding  experimental 
observation  (right),  for  the  thick-box,  multi-source  evaluation. 
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Appendix  B:  MCNP6  Input  Files 


This  appendix  presents  and  discusses  three  examples  of  the  30  MCNP6  input  decks  that 
produced  the  modeled  spectra  for  this  work.  Each  example  illustrates  key  elements  of  a 
particular  series  of  simulations. 

Initial  Validation 

c  cell  cards 

100  1  -0.001205  -10  20  30  40  50  60  imp:p=l  $  air  in  world  (rpp) 

200  2  -7.82  -20 : -3 0 : -40 : -5 0 : -60  70  imp:p=l  $  iron  shield  (lots  of  rpp's)  1  -0.001205 

300  1  -0.001205  -70  80  imp:p=l  $  air  in  shielded  box 

400  3  -3.67  -80  imp:p=l  $  Nal  detector 

999  0  10  imp:p=0  $  rest  of  the  world 


c  surface  cards 


10  rpp 
20  rpp 
30  rpp 
40  rpp 
50  rpp 
60  rpp 
70  rpp 
80  rcc 


0  255.88  0  255.88  0  255.88 


99.365 

156.515  99.365  156.515  99.365  156.515 

97.46 

158.42 

97.46 

158.42 

155.88 

156.515 

97.46 

158.42 

97.46 

158.42 

99.365 

100 

97.46 

158.42 

99.365 

100 

99.365 

156.515 

97.46 

158.42 

155.88 

156.515  99.365  156.515 

100 

155.88 

100  155.88 

100  155.88 

127.94 

127.94 

124.13  0  0  7.62 

3.81 

$  problem  bounding  surface 
$  steel  box  cube 
$  steel  box  lid 
$  steel  box  base 
$  steel  box  left  side 
$  steel  box  right  side 
$  air  inside  the  box 
$  Nal  Detector,  cylinder 


c  material  specification 
ml  07014  0.8  8016  0.2 
m2  26000  0.977170  6012  0.022831 
m3  11023  -0.5  53127  -0.5 
mode  p 

SDEF  POS=90  127.94  127.94  PAR=p  ERG=.662 
.6617  MeV 
nps  1000000 
f8:p  400 

E8  0  l.e-5  .01  .02  .03  .04  .05  .06  .07  .08  .09  & 

.1  .11  .12  .13  .14  .15  .16  .17  .18  .19  .20  .21  .22  . 
.26  .27  .28  .29  & 

.30  .31  .32  .33  .34  .35  .36  .37  .38  .39  & 

.40  .41  .42  .43  .44  .45  .46  .47  .48  .49  & 

.50  .51  .52  .53  .54  .55  .56  .57  .58  .59  & 

.60  .61  .62  .63  .64  .65  .66  .662 


$  air 

$  shield,  carbon  steel 
$  Sodium  Iodide 

$  source  position,  type,  energy  photons  of  energy 

$  number  of  photons  generated 
$  Energy  deposition  in  detector 
$  Energy  bins  for  the  spectrum  in  1  Oke V  increments 
.24  .25  & 


The  shielded  and  unshielded  input  decks  are  nearly  identical,  and  this  remained  true  for 
all  subsequent  simulations  that  produced  a  GPF.  Swapping  the  highlighted  portions 
changed  the  “box  cell”  material  and  density  from  air  to  iron.  This  technique  is  effectively 
an  on/off  switch  for  the  simulation’s  shield  for  every  simulated  GPF  produced. 


Alternate  Source  Study 

c  cell  cards 

100  1  -0.001205  -10  20  imp:p=l  $  air  in  world  (rpp) 
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200 

1 

-0.001205  -20  70 

imp:p=l 

$  iron  shield  (lots  of  rpp's)  2  -7.82 

300 

1 

-0.001205  -70  80 

imp:p=l 

$  air  in  shielded  box 

400 

3 

-3.67  -80 

imp:p=l 

$  HPGe  detector 

999 

0 

10 

imp:p=0 

$  rest  of  the  world 

c  surface  cards 

10  rpp  0  255.88  0  255.88  0  255.88 

20  rpp  99.365  156.515  99.365  156.515  99.365  156.515 
70  rpp  100  155.88  100  155.88  100  155.88 

80  rcc  127.94  127.94  123.525  0  0  8.83  3.15 

c  material  specification 
ml  07014  0.8  8016  0.2  $  air 

m2  26000  0.977170  6012  0.022831  $  shield,  carbon  steel 

m3  11023  -0.5  53127  -0.5  $  Sodium  Iodide 

mode  p 

SDEF  POS=92.38  127.94  127.94  PAR=p  ERG=D1  $  source  position,  type,  energy 
SI1  L  1.173  1.332 
SP1  1.0  1.0 

nps  10000000  $  number  of  photons  generated 

f8:p  400  $  Energy  deposition  in  detector 

e8  0  l.e-5  & 


0 

.01  0 

.02  0. 

,03  0 

.04  0 

.05  0 

.06  0 

.07  0 

.08  0 

.09  & 

0.10 

0.11 

0.12 

0.13 

0.14 

0.15 

0.16 

0.17 

0.18 

0.19 

& 

0.20 

0.21 

0.22 

0.23 

0.24 

0.25 

0.26 

0.27 

0.28 

0.29 

& 

0.30 

0.31 

0.32 

0.33 

0.34 

0.35 

0.36 

0.37 

0.38 

0.39 

& 

0.40 

0.41 

0.42 

0.43 

0.44 

0.45 

0.46 

0.47 

0.48 

0.49 

& 

0.50 

0.51 

0.52 

0.53 

0.54 

0.55 

0.56 

0.57 

0.58 

0.59 

& 

0.60 

0.61 

0.62 

0.63 

0.64 

0.65 

0.66 

0.67 

0.68 

0.69 

& 

0.70 

0.71 

0.72 

0.73 

0.74 

0.75 

0.76 

0.77 

0.78 

0.79 

& 

0.80 

0.81 

0.82 

0.83 

0.84 

0.85 

0.86 

0.87 

0.88 

0.89 

& 

0.90 

0.91 

0.92 

0.93 

0.94 

0.95 

0.96 

0.97 

0.98 

0.99 

& 

1.00 

1.01 

1.02 

1.03 

1.04 

1.05 

1.06 

1.07 

1.08 

1.09 

& 

1.10 

1.11 

1.12 

1.13 

1.14 

1.15 

1.16 

1.17 

1.18 

1.19 

& 

1.20 

1.21 

1.22 

1.23 

1.24 

1.25 

1.26 

1.27 

1.28 

1.29 

& 

1.30 

1.31 

1.32 

1.33 

1.34 

This  series  of  simulations  varied  the  sources  to  simulate  the  emission  spectra  of  60Co, 
57Co,  and  109Cd.  All  sources  were  simulated  as  point  sources.  The  highlighted  portion  is 
the  section  varied  to  produce  the  appropriate  source  emissions,  which  in  this  case 
matches  the  energies  and  relative  activities  of  60Co. 


$  problem  bounding  surface 
$  steel  box  cube 
$  air  inside  the  box 
$  Nal  Detector,  cylinder 


Thick  Wall  Assembly 

c  cell  cards 

100  1  -0.001205  -30  20 

200  2  -7.874  -20  70 

300  1  -0.001205  -70  80 

400  3  -3.67  -80 

500  4  -2.4  -10  30 

999  0  10 


imp:p=l 

imp:p=l 

imp:p=l 

imp:p=l 

imp:p=l 

imp:p=0 


$  air  in  room  (rpp) 

$  iron  shield  (rpp)  2  -7.874 

$  air  in  shielded  box 
$  Nal  detector 

$  Concrete  floor,  walls,  and  ceiling 
$  rest  of  the  world 


c  surface  cards 
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10 

20 

30 

70 

80 


rpp 

rpp 

rpp 

rpp 

rcc 


-100  355.88  -100  355.88  -100  355.88 
96.825  159.055  96.825  159.055  96.825  159. 
0  255.88  0  255.88  0  255.88 

100  155.88  100  155.88  100  155.88 

127.94  127.94  124.13  0  0  7.62  3.81 


$  problem  bounding  surface 
055  $  Outer  Steel  Surface 
$  Inner  surface  of  concrete  enclosure 
$  Inner  steel  surface 
$  Nal  Detector,  cylinder 


c  material  specification 
ml  07014  0.778  07015  0.00286  & 

8016  0.20896  8017  0.000817  8018  0.000421  & 

18040  0.00930  18036  0.000315  18038  0.00005884 
m2  26056  0.896  26054  0.0567  & 

26057  0.0215  26058  0.00274  & 

6012  0.0226  6013  0.000251 
m3  11023  -0.5  53127  -0.5 
m4  1001  -0.01  6012-0.001  8016  -.529107  & 

11023  -0.002  12024  -0.002  13027  -0.033872  & 

14028  -0.337021  19039  -0.013  20040  -0.044  & 

26000  -0.14000  SConcrete 

mode  p 

sdef  x=92.38  y=92.38  z=D20 

erg=FZ=D45  $Position  and  energy 

si20  L  127.94  163.5  163.5  $  Discrete  Z  coordinates  for  point  sources 

sp20  254.8685  147.7268  147.7268  $  Relative  Probability 

ds45  L  0.662  1.173  1.332  $  Energies  in  MeV 

nps  10000000  $  number  of  photons  generated 

f8:p  400  $  Energy  deposition  in  detector 

e8  0  l.e-5  & 


0 

.01  0 

.02  0. 

,03  0 

.04  0 

.05  0 

.06  0 

.07  0, 

.08  0 

.09  & 

0.10 

0.11 

0.12 

0.13 

0.14 

0.15 

0.16 

0.17 

0.18 

0.19 

& 

0.20 

0.21 

0.22 

0.23 

0.24 

0.25 

0.26 

0.27 

0.28 

0.29 

& 

0.30 

0.31 

0.32 

0.33 

0.34 

0.35 

0.36 

0.37 

0.38 

0.39 

& 

0.40 

0.41 

0.42 

0.43 

0.44 

0.45 

0.46 

0.47 

0.48 

0.49 

& 

0.50 

0.51 

0.52 

0.53 

0.54 

0.55 

0.56 

0.57 

0.58 

0.59 

& 

0.60 

0.61 

0.62 

0.63 

0.64 

0.65 

0.66 

0.67 

0.68 

0.69 

& 

0.70 

0.71 

0.72 

0.73 

0.74 

0.75 

0.76 

0.77 

0.78 

0.79 

& 

0.80 

0.81 

0.82 

0.83 

0.84 

0.85 

0.86 

0.87 

0.88 

0.89 

& 

0.90 

0.91 

0.92 

0.93 

0.94 

0.95 

0.96 

0.97 

0.98 

0.99 

& 

1.00 

1.01 

1.02 

1.03 

1.04 

1.05 

1.06 

1.07 

1.08 

1.09 

& 

1.10 

1.11 

1.12 

1.13 

1.14 

1.15 

1.16 

1.17 

1.18 

1.19 

& 

1.20 

1.21 

1.22 

1.23 

1.24 

1.25 

1.26 

1.27 

1.28 

1.29 

& 

1.30 

1.31 

1.32 

1.33 

1.34 

$  air 


$  shield,  carbon  steel 
$  Sodium  Iodide 


This  series  of  simulations  increased  the  thickness  of  the  steel  enclosure  to  1.25”.  It  also 
defines  a  small  concrete  room  around  the  box.  The  location  of  the  source  was  then  varied, 
and  the  type  of  source  was  interchanged  between  Cs  and  C0-60.  This  particular 
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example  is  the  multi-source  configuration,  which  has  both  types  of  sources:  Cs  outside 
an  edge,  and  60Co  outside  a  comer. 


Fidelity  Study 

The  example  above  also  shows  some  of  the  variations  that  came  from  the  fidelity  study. 
This  series  of  simulations  varied  the  material  data  cards,  the  number  of  particles  in  the 
source  card,  and  the  surface  and  cell  cards  defining  the  geometry.  One  noteworthy 
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difference  in  the  example  provided  is  the  complexity  of  the  material  data  cards,  which 
include  more  elements  and  isotopes  than  the  earlier  series  of  models. 
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Appendix  C:  MATLAB  Monte  Carlo  Script 


%This  is  the  Monte  Carlo  code  for  transporting  monoenergetic  photons  from 
%a  point  source  to  a  detector  between  two  iron  barrier. 

clear  all 
clc 

%This  simulation  tracks  photons  as  they  Compton  scatter  through  a  slab 
“/..geometry  configuration  to  represent  a  point  source  . 

“/..INPUTS 

%To  change  the  setup,  you  should  alter  the  following  variables 
tracks=10000000;  %Number  of  primary  photon  tracks  to  be  evaluated. 

primeenergy=.662;  %Primary  photon  energy  in  MeV,  from  .001  to  20 


% - 

%SURFACE  LOCATION  INPUTS 

%Below  is  a  crude  diagram  of  the  setup  showing  where  each  surface  boundary 
%is  located.  These  boundaries  can  be  changed,  but  be  careful! 

%Each  boundary  must  be  given  a  coordinate  relative  to  the  other  boundaries 
%A11  coordinates  are  in  cm. 


%  negair 

src  nio 

nii  Nalpos  fli 

% 

%  | 

s 

|  iron  | 

Nal 

|  iron  | 

%  | 

s 

|  iron  | 

Nal 

|  iron  [ 

%  | 

s 

|  iron  | 

Nal 

|  iron  | 

%  |  air 

s 

air  |  iron 

air  Nal 

air  |  iron  | 

%  | 

s 

|  iron  | 

Nal 

|  iron  | 

%  | 

s 

|  iron  | 

Nal 

|  iron  | 

%  | 

s 

|  iron  | 

Nal 

|  iron  | 

src  =  0;  “/..position  of  the  source 

negair=  -100;  %left  boundary  of  the  simulation,  must  be  negative, 
nio  =  10;  “/..outer  surface  of  the  near  iron  barrier,  must  be  positive, 
nii  =  11;  “/..inner  surface  of  the  near  iron  barrier,  must  be  greater  than  nio. 
Nalpos  =  41.5;  “/..position  of  center  of  the  Nal  detector, 
fli  =  70;  “/..inner  surface  of  the  far  iron  barrier,  must  be  greater  than  fNal. 
fio  =  71;  “/..outer  surface  of  the  far  iron  barrier  and  far  boundary  of  the 
“/..simulation,  must  be  the  highest  of  these  coords. 

% - 


%We  begin  by  defining  variables  that  will  endure  throughout  the  simulation. 

DetectorPosition=10000;  %This  shows  a  position  for  evaluating  the  flux. 

Energydata=zeros(l,  (tracks/2));  “/..Initialize  what  will  be  an  allocatable  array.  This  will 
“/..provide  the  photon  flux  output  of  the  Monte  Carlo  simulation. 
n=l;  %This  is  a  counter  for  the  number  of  photons  that  pass  through  the 
“/..detector's  coordinate. 

% 

%Now  we  input  the  NIST  data  tables 

AirEnergies  =  [1.00000E-03,  1.50000E-03,2.00000E-03,3.00000E-03,  3.20290E-03,  3.202901E- 

03,4.00000E-03,5.00000E-03,6.00000E-03,8.00000E-03,1.00000E-02,1.50000E-02,2.00000E-02,3.00000E 

02,4.00000E-02,5.00000E-02,6.00000E-02,8.00000E-02,1.00000E-01,1.50000E-01,2.00000E-01,3.00000E 
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01,4.00000E-01,5.00000E-01,6.00000E-01,8.00000E- 

01, 1.00000E+00,1.25000E+00,1.50000E+00,2.00000E+00,3.00000E+00,4.00000E+00,5.00000E+00, 6.00000 
E+00,8.00000E+00,1.00000E+01,1.50000E+01,2.00000E+01]; 

AirCoefficients= 

[3.606E+03,1.191E+03,5.279E+02,1.625E+02,1.340E+02,1.485E+02,7.788E+01,4.027E+01,2.341E+01,9.9 

21E+00,5.120E+00,1.614E+00,7.779E-01,3.538E-01,2.485E-01,2.080E-01,1.875E-01,1.662E-01,1.541E- 

01,1.356E-01,1.233E-01,1.067E-01,9.549E-02,8.712E-02,8.055E-02,7.074E-02,6.358E-02,5.687E- 

02,5.175E-02,4.447E-02,3.581E-02,3.079E-02,2.751E-02,2.522E-02,2.225E-02,2.045E-02,1.810E- 

02,1.705E-02]; 

IronEnergies  =  [.001, 1.5*10A-3,  2*10A-3,  3*10A-3,  4*10A-3,  5*10A-3,  6*10A-3,  7.111*10A-3,  7.112*10A-3, 
8*10A-3, 10*10A-3, 15*10A-3,  20*10A-3,  30*10A-3,  40*10A-3,  50*10A-3,  60*10A-3,  80*10A-3, 100*10A-3, 
150*10A-3,  200*10A-3,  300*10A-3,  400*10A-3,  500*10A-3,  600*10A-3,  800*10A-3, 1000*10A-3, 1250*10A- 
3, 1500*10A-3,  2000*10A-3,  3000*10A-3,  4000*10A-3,  5000*10A-3,  6000*10A-3,  8000*10A-3, 10000*10A-3, 
15000*10A-3,  20000*10A-3]; 

IronCoefficients  =  [9.085*10A2*10,  399*10A2*10, 1.626*10A2*10,  5.576*10A1*10,  2.567*10A1*10, 
1.398*10A1*10,  8.484*10,  5.319*10,  4.076*10A1*10,  3.056*10A1*10,  1.076*10A1*10,  5.708*10,  2.568*10, 
8.176*10A-1*10,  3.629*10A-1*10,  1.958*10A-1*10, 1.205*10A-1*10,  5.952*10A-2*10,  3.717*10A-2*10, 
1.964*10A-2*10,  1.460*10A-2*10,  1.099*10A-2*10,  9.4*10A-3*10,  8.414*10A-3*10,  7.704*10A-3*10, 
6.699*10A-3*10,  5.995*10A-3*10,  5.350*10A-3*10,  4.883*10A-3*10,  4.265*10A-3*10,  3.621*10A-3*10, 
3.312*10A-3*10,  3.146*10A-3*10,  3.057*10A-3*10,  2.991*10A-3*10,  2.994*10A-3*10,  3.092*10A-3*10, 
3.224*10A-3*10]; 


denair=.0012922;  %density  of  air,  in  grams/cmA3 
deniron  =  7.874;  %density  of  iron,  in  grams/cmA3 

%The  Monte  Carlo  program  is  governed  by  a  for  loop  that  starts  a  primary 
%photon  and  runs  it  through  a  set  number  of  tracks.  These  tracks  continue 
%through  multiple  interactions  until  a  kill  criteria  is  met,  at  which  point 
%a  new  primary  photon  is 
%created. 
for  i  =  l:tracks 

location  =  src;  %This  is  the  source  location,  which  can  be  moved,  in  cm. 
energy  =  primeenergy; 

angle  =  pi*rand(l);  %Initial  angle  of  the  photon  with  respect  to  the 
%plane  of  the  iron  barrier. 

while  (negair  <=  location  <=  fio)  %This  is  the  primary  loop  that  drives  each  track 
%It  will  keep  running  as  long  as 
%the  photon  remains 


mfp=(-log(l-rand(l)))/(log(2));  %Samples  the  number  of  mean  free  paths 
%a  particle  will  travel  before  interacting 


% - 

%This  section  produces  the  mean  free  path  length  for  air,  iron,  and  Nal 
%for  the  current  energy  of  the  photon.  This  uses  a  log-log  interpolation  of 
%NIST  data  for  air,  iron,  sodium,  and  iodine. 


%  We  start  with  two  lists  of  photon  energies  and  attenuation  constants 

%from  NIST. 

k=length(AirEnergies); 

while  k  >=  1  %This  finds  the  NIST  entry  just  below  the 

if  (AirEnergies(k)  >  energy)  %energy  of  the  photon. 

k=  k-1; 
else 
break 
end 


102 


end 

%Now  we  do  the  log-log  interpolation. 

lair  =  log(AirCoefficients(k))  *  (log(AirEnergies(k+l))  -  log(energy))  /  (log(AirEnergies(k+l))  - 
log(AirEnergies(k)))  +  log(AirCoefficients(k+l))  *  (log(energy)  -  log(AirEnergies(k)))  / 
(log(AirEnergies(k+l))  -  log(AirEnergies(k))); 

airAtten=exp(lair);  %This  is  the  cross  section  per  density,  cmA2  per  gram 
pathAir=(-log(0.5))/((airAtten)*denair); 

%Now  we  repeat  the  same  process  for  iron. 
k=length(IronEnergies); 

while  k  >=  1  %This  finds  the  NIST  entry  just  below  the 

if  (IronEnergies(k)  >  energy)  %energy  of  the  photon. 

k=  k-1; 
else 
break 
end 
end 

%Now  we  do  the  log-log  interpolation. 

liron  =  log(IronCoefficients(k))  *  (log(IronEnergies(k+l))  -  log(energy))  /  (log(IronEnergies(k+l))  - 
log(IronEnergies(k)))  +  log(IronCoefficients(k+l))  *  (log(energy)  -  log(IronEnergies(k)))  / 
(log(IronEnergies(k+l))  -  log(IronEnergies(k))); 

ironAtten=exp(liron);  %This  is  the  cross  section  per  density,  cmA2  per  gram 
pathIron=(-log(0.5))/((ironAtten)*deniron); 


% - 

%This  section  transports  the  photon  along  its  trajectory  until  its  number 
%of  mean  free  paths  has  been  reduced  to  zero. 


previouslocation=location;  %This  stows  the  starting  location  for  use  later. 

while  mfp  >  0  %This  loop  will  transport  the  photon  left  or  right  as  necessary  through  the  surfaces  of  the 
problem. 

if  location  <=  nio  %This  is  the  outer  air  loop.  It  identifies  if  the  photon  is  currently  in  the  air  outside 
the  box  and  transports  it. 

if  ((mfp*pathAir*cos(angle))+location)  >=  nio;  %If  it  moves  out  of  the  air  and  into  the  near  iron 
barrier 

mfp=mfp-(((nio-location)/(pathAir*cos(angle))));  %  reduces  mfp  appropriately 
location  =  nio;  %Now  it's  out  of  the  air  and  ready  for  the  near  iron  barrier  loop, 
else 

location  =  location  +  mfp*pathAir*cos(angle);  %shows  where  the  next  interaction  happens  if  it 
starts  and  ends  in  the  outer  air. 
mfp=0; 
end 
end 


if  nii>=location 

if  location  >=nio  %  These  if  statements  identify  if  the  photon  is  currently  inside  the  near  iron 
barrier 

if  ((mfp*pathIron*cos(angle))+location)  >=  nii;  %If  it  moves  into  the  air  inside  the  box 
mfp=mfp-(((nii-location)/(pathIron*cos(angle))));  %  reduces  mfp  appropriately 
location  =  nii;  %Now  it's  at  the  inner  iron  surface,  ready  for  the  inner  air  loop. 
elseif((mfp*pathIron*cos(angle))+location)  <=  nio;  %If  it  moves  into  the  air  outside  the  box 
mfp=mfp-abs((((location-nio)/(pathIron*cos(angle)))));  %  reduces  mfp  appropriately 
location  =  nio;  %Now  it's  outside  the  box,  ready  for  the  outer  air  loop, 
else 

location  =  location  +  mfp*pathIron*cos(angle);  %shows  where  the  next  interaction  happens  if 
the  beginning  and  end  are  in  the  air. 
mfp=0; 
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end 

end 

end 

if  fii>=location 

if  location  >=nii  %  These  if  statements  identify  if  the  photon  is  currently  in  the  'inner  air'  between 
the  iron  barriers 

if  ((mfp*pathAir*cos(angle))+location)  >=  fii;  %If  it  moves  into  the  far  barrier 
mfp=mfp-(((fii-location)/(pathAir*cos(angle))));  %  reduces  mfp  appropriately 
location  =  fii;  %Now  it's  at  the  inner  iron  surface,  ready  for  the  far  iron  loop. 
elseif((mfp*pathAir*cos(angle))+location)  <=  nii;  %If  it  moves  into  near  iron  barrier 
mfp=mfp-abs((((location-nii)/(pathAir*cos(angle)))));  %  reduces  mfp  appropriately 
location  =  nii;  %Now  it's  outside  the  box,  ready  for  the  outer  air  loop, 
else 

location  =  location  +  mfp*pathAir*cos(angle);  %shows  where  the  next  interaction  happens  if 
the  beginning  and  end  are  in  the  air. 
mfp=0; 
end 
end 
end 

if  fio>=location 

if  location  >=fii  %  These  if  statements  identify  if  the  photon  is  currently  inside  the  far  iron  barrier 
if  ((mfp*pathIron*cos(angle))+location)  <=  fii;  %If  it  moves  into  the  air  inside  the  box 
mfp=mfp-abs((((location-fli)/(pathIron*cos(angle)))));  %  reduces  mfp  appropriately 
location  =  fii;  %Now  it's  out  of  the  iron  and  ready  for  the  inner  air  loop. 

else 

location  =  location  +  mfp*pathIron*cos(angle);  %shows  where  the  next  interaction  happens  if 
the  beginning  and  end  are  in  the  iron. 
mfp=0; 
end 

end 

end 

end 

% - 

%The  photon  is  now  at  its  interaction  location.  Now  we  run  our  kill  criteria. 

if  (location)  <  negair  %This  is  a  kill  criteria,  any  photon  that  makes  it  a  meter  away  from  the  source  in 
the  wrong  direction  is  killed. 

break 

end 

if  (location)  >  fio  %This  is  a  kill  criteria,  any  photon  that  makes  it  through  the  far  side  iron  barrier  is 
killed. 

break 

end 

% - 

%Now  we  must  determine  if  the  photon  passed  through  the  plane  of  interest. 

if  previouslocation  <=  Nalpos  %Checks  to  see  if  the  photon  has  crossed  the  altitude  of  interest  from  left 
to  right 

if  location  >  Nalpos 

Energydata(n)=energy;  %adds  the  photon  energyto  the  altitude  of  interest. 
n=n+l; 
end 
end 

if  previouslocation  >=  Nalpos  %Checks  to  see  if  the  photon  has  crossed  the  altitude  of  interest  from 
right  to  left 
if  location  <  Nalpos 
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Energydata(n)=energy;  %adds  the  photon  energyto  the  altitude  of  interest.  %adds  the  photon 
energy 

%to  the  altitude  of  interest. 
n=n+l; 
end 
end 

% - 

%  The  following  section  establishes  and  samples  from 

%the  Klein  Nishina  distribution.  I  use  rejection  sampling  for  this. 

sample  =  0; 

gamma=energy/.5 1 1 ; 

while  sample  <  1 

fcompmax  =  energy/(l+gamma*(l-cos(0))); 

knprobmax  =  (l/2)*((fcompmax)A2)*(l+(cos(0))A2+((((gamma)A2)*((l-cos(0))A2))/(l+gamma*(l- 
cos(0))))); 

%(l/2)*((fcompmax)A2)*(fcompmax+((fcompmax)A(-l))-((sin(fcompmax))A2));  %Normalizes  the 
probability  distribution. 
knsample=pi*rand(l); 

fcomp  =  energy/(l+(gamma)*(l-cos(knsample))); 
knprob  =  (l/2)*((fcomp)A2)*(l+(cos(knsample))A2+((((gamma)A2)*((l- 
cos(knsample))A2))/(l+gamma*(l-cos(knsample))))); 

%(l/2)*((fcomp)A2)*(fcomp+((fcomp)A(-l))-((sin(knsample))A2)); 

if  rand(l)<=(knprob/knprobmax)  % 
scatterangle=knsample; 

sample=2;  %The  sample  is  accepted,  this  value  kills  the  while  loop, 
end 
end 


if  energy  <  0.1  %Below  lOOkeV  in  air,  the  photoelectric  effect  becomes 
%dominant  in  iron.  This  loop  kills  low-energy 
%photons  that  interact  in  iron. 

if  location  <=  nio  %This  loop  kills  the  track  for  low-energy 
“/□interactions  inside  the  near  iron  barrier, 
if  location  >=  nii 
break 
end 
end 

if  location  <=  fli  %This  loop  kills  the  track  for  low-energy 
“/□interactions  inside  the  far  iron  barrier, 
if  location  >=  fio 
break 
end 
end 
end 


%Now  we  need  the  energy  of  the  scattered  photon.  The  Compton  Formula 
%provides  this  from  the  scattered  angle. 
energy=energy/(l+energy*(l-cos(scatterangle))); 

if  energy  <.005  %This  loop  is  an  energy-based  kill  criteria  for  the  track. 

break 

end 


%■ 
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%We  cannot  simply  add  the  scattering  angle  to  the  angle  of  the  photon. 

%This  is  because  a  scattering  angle  is  equally  likely  all  around  the  axis 
%of  the  incident  photon.  Thus,  while  we  could  simply  add  it  to  the  zero 
%angle  of  the  initial  photon,  we  cannot  do  so  for  further  interactions. 

%Since  we  are  only  interested  in  the  projection  of  the  scattered  photon 
%onto  the  axis  of  the  plane  of  the  slab,  we  can  use  the  Law  of  Cosines  for 
%Spherical  Triangles. 

revolvingangle  =  2*pi*rand(l); 

angle  =  acos(cos(angle)*cos(scatterangle)+sin(angle)*sin(scatterangle)*cos(revolvingangle)); 

end 

end 

Spectrumdata=Energydata(Energydata~=0); 

%Dosedata=zeros(l,  (length(Spectrumdata)));  %Initialize  what  will  be  an  allocatable 
%array.  This  will  track  dose  deposited  in  the  detector  volume. 

%for  m=l:length(Spectrumdata) 

%  Dosedata(m)=(Spectrumdata(m)*(1.6*10A(-13))*NaAtten)+(Spectrumdata(m)*(1.6*10A(- 
13))*IAtten); 

%end 

%Dose=sum(Dosedata) 

hist(Spectrumdata,662) 
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